{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true,
    "slideshow": {
     "slide_type": "skip"
    }
   },
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "%config InlineBackend.figure_format='png'\n",
    "import pylab as pl\n",
    "import numpy as np\n",
    "import scipy.signal as signal\n",
    "import matplotlib as mpl\n",
    "mpl.rcParams['axes.grid'] = True\n",
    "import matplotlib.ticker as ptick "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "## Weaver方式SSB復調のためのIIRフィルタ設計"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "NanoSDR(FriskSDR)のため設計したフィルタを使用したが阻止量(40dB)に不足を感じたので再設計する。\n",
    "fs=48kHzで、fc=1300Hzは同じ。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "### 方針と仕様\n",
    "\n",
    "2.6kHz帯域のSSBをWeaver方式で復調するため、0Hzを中心とした正負対称で1300HzカットオフのLPFを設計する。\n",
    "\n",
    "* サンプリングレート48000Hzにてfc=1300Hzの6次IIRローパスフィルタ\n",
    "* 楕円(ellipitic)フィルタで、リプル1dB, 阻止量80dBとする\n",
    "* scipy.signalに用意されているフィルタ設計関数を用いて設計する"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([ 0.00991993, -0.05629584,  0.13612359, -0.17949116,  0.13612359,\n",
       "        -0.05629584,  0.00991993]),\n",
       " array([  1.        ,  -5.79238023,  14.03316532, -18.19974154,\n",
       "         13.32534556,  -5.22213803,   0.85575361]))"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAD7CAYAAAB5aaOHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAHzBJREFUeJzt3XmQHOWZ5/Hvo6Nl1AK1JIPAakMLi0PChgaDwCdtGGTh\nCI5lbBnPGNxj8GCDwYyx12AToTkIYw94EcbLTsQurIbDKGyHB8vHgoRRwSwhrgVxSYbmEKDWARgk\nsEBX69k/3iq9pdZV3VVvV1Xm7xNRUfVmZmW/+ahUT71HZpq7IyIisifD6l0BERFpDkoYIiJSESUM\nERGpiBKGiIhURAlDREQqooQhIiIVGVHvCuyOmWnOr4jIILi71XqfDd/CcHc93Jk9e3bd69AoD8VC\nsVAsdv9IpeEThgTLly+vdxUahmIRKRaRYpGeEoaIiFRECaNJdHd317sKDUOxiBSLSLFIz1L2d1XL\nzLyR6yci0ojMDM/joLcEhUKh3lVoGIpFpFhEikV6ShgiIlIRdUmJiGSMuqRERKSulDCahPpnI8Ui\nUiwixSI9JQwREamIxjBERDJGYxgiIlJXShhNQv2zkWIRKRaRYpFe3RKGmc00sz+Z2XNm9r161UNE\nRCpTlzEMMxsGPAecDKwEHgHOdvc/9dtOYxgiIgOUtTGM6UCPu7/s7puBecAZdaqLiIhUoF4JYxLw\nall5RXGZ7IL6ZyPFIlIsIsUivYa+RWt/S5bAvHmwfj18+tNw1lkwfPj22zz6KCxYABdeCG1t9amn\niEgW1Sth9AIHlpXbi8t2cO653WzY0MFDD8HatW2ccUYnRx/dxZw5cNllBS64AL7//S7M4KqrClxz\nDZxwQhdLlsCddxa47TaYNasLiL9Aurqar9zV1dVQ9VG5ccoljVKfepVLyxqlPkNZLhQKzJ07F4CO\njg5Sqdeg93DgWcKg9yrgYeBL7r6s33Y+YYIzbRpcfDGceSaMHBnWucPdd8Nll8HYsXDAAbB4Mcyf\nD0ccAcccA3/6EyxaBGWfJxGRzMvUoLe79wHfBBYAzwDz+ieLksceg/vvhy98ISYLADOYOROeeAKu\nuCIkk2eegWOPhb32gj/8ISSNlSuH4ojS6/9rMs8Ui0ixiBSL9Oo2huHudwGH7Wm7Aw/c/foRI+C0\n03ZcPnkynHhidhKGiEi9ZfpaUtdeC729cN11NayUiEiDy1SX1FCZNEktDBGRWsl0wjjggOwkDPXP\nRopFpFhEikV6mU4Yra3w3nv1roWISDZkegzjiSfgnHPgySdrWCkRkQanMYxBaGmBzZvrXQsRkWzI\nfMLYtKnetagN9c9GikWkWESKRXqZThgjR2YnYYiI1FumxzBWr4ajjoI1a2pYKRGRBqcxjEHIUpeU\niEi9ZT5hZGXQW/2zkWIRKRaRYpFe5hOGWhgiIrWR6TEMdxg2DLZuDVe3FRHJA41hDIJZuJptVrql\nRETqKdMJA7LTLaX+2UixiBSLSLFITwlDREQqkukxDICJE8O1pCZOrFGlREQanMYwBkktDBGR2lDC\naBLqn40Ui0ixiBSL9DKfMHQ9KRGR2sj8GMaRR8Ktt4ZrSomI5IHGMAYpK11SIiL1louEkYUT99Q/\nGykWkWIRKRbp5SJhqIUhIlK9zI9hzJgB3/lOeBYRyQONYQySZkmJiNRG5hNGVrqk1D8bKRaRYhEp\nFunlImFkYdBbRKTeMj+Gcc45YfzinHNqVCkRkQanMYxBykqXlIhIvSlhNAn1z0aKRaRYRIpFeplP\nGJolJSJSG8nGMMzsX4HTgI3AC8DfufvbxXVXAF8FtgDfcvcFu9hH1WMY3/lOuBfGd79b1W5ERJpG\nM45hLACOcPdOoAe4AsDMpgGzgKnAqcCNZlbzAyvRLCkRkdpIljDc/R5331osPgi0F1+fDsxz9y3u\nvpyQTKanqofGMLJHsYgUi0ixSG+oxjC+Cvyh+HoS8GrZut7isiSykjBEROptRDVvNrOFQPndsg1w\n4Afu/tviNj8ANrv7HYP5G93d3XR0dADQ1tZGZ2cnXV1dQPxFsbvyq69Ca2vl2zdquaurq6Hqo3Lj\nlEsapT71KpeWNUp9hrJcKBSYO3cuwLbvyxSSnrhnZt3A14CT3H1jcdnlgLv7j4vlu4DZ7v7QTt5f\n9aD39dfDCy/AT39a1W5ERJpG0w16m9lM4LvA6aVkUTQfONvMWsxsMjAFeDhVPVoy0iXV/9dknikW\nkWIRKRbpVdUltQc3AC3AwuIkqAfd/UJ3X2pmvwCWApuBC6tuRuxGi2ZJiYjUROavJXXrrbBgQXgW\nEcmDpuuSahRZ6ZISEak3JYwmof7ZSLGIFItIsUgv8wlD15ISEamNzI9hLFgA11wDCxfWqFIiIg1O\nYxiDpFlSIiK1kYuEkYUuKfXPRopFpFhEikV6ShgiIlKRzI9hPPkkfPnL4VlEJA80hjFImiUlIlIb\nmU8YWRn0Vv9spFhEikWkWKSXi4ShFoaISPUyP4axZg0ceWR4FhHJA41hDJJaGCIitZH5hJGVQW/1\nz0aKRaRYRIpFeplPGGphiIjURubHMNxh2DDYuhWs5j16IiKNR2MYg2SWnW4pEZF6ynzCABgxAvr6\n6l2L6qh/NlIsIsUiUizSy0XCGD4ctmypdy1ERJpb5scwAMaNgxdfDM8iIlmnMYwqqIUhIlK9XCQM\njWFki2IRKRaRYpFeLhKGWhgiItXLxRjGgQfCf/4nHHRQDSolItLgNIZRhREj1MIQEalWLhLG8OEa\nw8gSxSJSLCLFIr1cJIwsDHqLiNRbLsYwPvxhuOMO+MhHalApEZEGpzGMKqiFISJSvVwkjCxMq1X/\nbKRYRIpFpFikl4uEoRaGiEj1cjGG8fGPwzXXwCc+UYNKiYg0uKYdwzCzy8xsq5mNL1t2hZn1mNky\nM5uRug5ZmFYrIlJvSROGmbUDpwAvly2bCswCpgKnAjeapb0XXhZO3FP/bKRYRIpFpFikl7qFcR3w\n3X7LzgDmufsWd18O9ADTU1ZCLQwRkeqNSLVjMzsdeNXdn+rXgJgELC4r9xaXJVPJoPfTT4frTX3o\nQ/BXfxXuA95Iurq66l2FhqFYRIOJhXu4x31fX3js6nWpvHVrfI/79q9rvay6/XRx223hdek4y495\nZ8+VLhvo9rXYx7Bh8K1v0VCqShhmthCYWL4IcOBK4PuE7qiqdHd309HRAUBbWxudnZ3b/pOUmqB7\nKg8f3sXPfgaLFxf4zGd2XP/0011cdRUce2yBnp7wwbvlFtiwobL9q5yN8r33Fti8GY4/vouNG8P6\nzZuhszOUFy8O5alTQ/nxxwts2QKHHNLFli3wzDMF+vqgo6OLzZvhuefC+g9+MJRfeimUJ04M27/y\nSth+woRQXrUqrB87Nmz/xhth/ZgxXfT1wbp1BbZuhb326mLrVnjnnVBuaQnld98N5eHDQ3njxvD+\nYcPC+zdvDuuhCzMwKzB8OIwc2cWwYeBeYNgwGDWqqzgVvYAZjB4d1m/YEMqtreH9pf8fY8aE9evX\nh/fvvXdYv3592H7s2FB+551QbmsL269bF8rjx4f1a9eG8oQJYf2bb4byvvuG9eXlYcNCfMxCPM3g\n9ddDef/9Q3n16lgGWLMm1PeAA0J59epYNgvxN9tx/Qc+EMqrVsWyGaxcGbYvrV+5Mqxvbw/l3t6w\nftKksH1v747rIXw+AFasCNu3t8ftjzqqss9voVBg7ty5ANu+L1NIMkvKzD4M3AO8S0gi7YSWxHTg\nqwDu/qPitncBs939oZ3spyazpM44A+bPD6/77+7xx+Gzn4VHHglXs3WH3/wGLrgArr0Wzjmn6j9f\nE4VCIbe/rN1h/Xp4++3wWLSowOGHd/Huu2H5u+/u/LGrdRs2hMfGjds/tmyBlhYYNaryx8iR8TFi\nxM6fB7qu/PXw4ds/hg3b/vWjjxb4+Me7dli3s23LX6cdNayPPP8f6S/VLKkkXVLu/jSwf6lsZi8B\nx7j7W2Y2H7jdzP4boStqCvBwinqUDB8eX7/zDuy9dyz/8Idw5ZXx0udmcOaZcNhhcMop4T/X3/5t\nytrlw6ZN8Oc/b/944434et26mBBKj9Kyd96B970P9tkHxo4NCeQDH4DRo3d8tLZCW9uO61tbw/Ne\ne4V97ezLv6Wl+b5IV62CyZPrXQvJiyE5D8PMXgSOdfc3i+UrgPOAzcC33H3BLt5XkxbGrFnwy1+G\n18uWweGHh9dvvAFTpsArr4Qvo/6WLoWurtDi+NjHqq5G5vT1wWuvhS+tlSu3fy49Xn89JIT33oPx\n42HChPh4//vj67a28G9Q/hg7NjzvvXf4xS0ilWmqFkZ/7n5wv/LVwNVD8bchtjD23ht6e2PCuOce\nOPHEnScLgGnTYO5c+Pzn4cEH4YMfHJLqNgx3WLMGXngBXnoJli8Pz6VHby+MGxd+zR9wQHw+8kiY\nOTO83nffkBD22af5fr2LyPZy8butlDAmTw5fciULFsCMPZw2+LnPwYUXQnc3LFxYv9lTKftn3WHF\nCnjqqdCqWrYsPJYuDf3oH/pQiN3kyXDCCfClL0FHR7iTYUtLkirtlvqqI8UiUizSy0XCKHVnTJ4c\nvhhLHn0ULrpoz++//HL4/e/hZz+DSy5JU8eh9Npr8NBD4fhLD4CjjgqtqunTQ4KcOjV0G4mIQE6u\nJXXeeXDzzXDxxaGFMGcObN4cuknefDMMhO5JT08Yx1i8GA45pOoqDak1a+C++6BQCI9Vq+D44+HY\nY8Pjox+F9nZ1GYlkRVOPYdRbqYXx/veHfniAZ58NM6MqSRYQksTll8Oll4bWRiNzh2eeCYP18+fD\nc8/Bpz4VBvDPPz+0JMpnjomIVKLBzmdOo/TluO++YZomhD76adMGtp9LLoHnn69PwqjkOjkvvQT/\n/M9w6KFw2mmh6+mHPwwtjPnz4dvfhmOOaf5koWsGRYpFpFikl4uEUd7CKCWMFSsGPuuppQWuvz60\nMjZtqm0dB6uvD+68M7Qepk8P01h//nN48cVQ15NPrs/AtIhkTy66pEq/qMsTRm8vTBrEFaxmzgyz\nhm6+Gb7+9drVcU/6z/7YsAFuugmuuy5MW/2Hf4CzzspHctBMmEixiBSL9HLbwujtDecNDMY//VPo\n6tm4sTb1G4i+vnBuyGGHwd13wy23hHNEzj47H8lCROonFwmjZNy46lsYEGYYffjDoZUxVAqFAo89\nFrqdbropdDvNnx/uJpi32U3qq44Ui0ixSC8XCaN0afN99gnXJ4JwCYvBJgyAf/zH0MoYirGMjRvh\n3/4NTj01DLzff79uNysiQy8XCSNc0hnGjAlXK+3rg9WrYf/9d/++3Zk+PUy1/Y//qE0dd2X58jAl\ndsOGLp56Cr7ylfy1KPpTX3WkWESKRXq5SBilFsawYeGqpevWhcQxZkx1+/3mN8PZ36ksXhwuxXH2\n2SEx7bdfur8lIrInuUgYpRYGhG6plStD4qj2ulCnnx5aAEuWVLefnfntb8N9PObODedP3HdfofZ/\npEmprzpSLCLFIr1cJIzy27OWEkb5PTEGa8QI+MY3at/KuOuucEb2738fpvGKiDSCXCSM8hZGS0u4\nP0MtEgaEL/Zf/Src5KcWHngAzj03nIx33HFxufpnI8UiUiwixSK9XCSM8hbGiBHw1lu1Sxj77RcG\npe+8s/p9rVwZ7r1xyy26YZOINJ7cJYyRI8MVamuVMCDcwvX226vbx5YtYXD7oot23g2l/tlIsYgU\ni0ixSE8JowZOPz3cX2LNmsHv46qrwpVzv//92tVLRKSWcnE/jL/5G7jjjnDZ75NPDneK27Sp+lZB\nuS9+Mdy977zzBv7epUvDrWKfeGLwlysRESlJdT+MXLQwyge9az2GUXLaaWEq7EC5h26o2bOVLESk\nseUiYaTukoJw2Y5774X33hvY+xYsCHfA29OVb9U/GykWkWIRKRbp5SJhlLcwUiWMCRPg6KND0qiU\nO1x5JfzLv8Qr6oqINKpcJIyTTgr3sICYMKq9LMjOzJgBf/xj5dsvWhQuUfLXf73nbTXHPFIsIsUi\nUizSy0XCuOiicGtVCL/k33kHRo+u/d858US4777Kt7/hBrj44uovUSIiMhRy91U1ciSsX5/mZkPH\nHQfPPhsvob47L78cLlP+5S9Xtm/1z0aKRaRYRIpFerlMGO4walTt9z1qVLjs+QMP7Hnb228PU3FT\ndI2JiKSQy4QB6W5nWmm31C9/CbNmVb5f9c9GikWkWESKRXq5Sxil2UipEsbHPgYPP7z7bZ5/Pkyl\n/dSn0tRBRCSF3CWMUgsjRZcUwFFHhftj7O4E9V/9Cs46C4YPr3y/6p+NFItIsYgUi/RymzBStTAm\nToT3vQ9efXXX29xzD3zuc2n+vohIKrlNGKlaGACdnbu+C9+mTeFChZ/85MD2qf7ZSLGIFItIsUgv\nacIws4vNbJmZPWVmPypbfoWZ9RTXzUhZh/5Sj2FA6JZ64omdr3vsMZgyBdra0v19EZEUkiUMM+sC\nTgM+4u4fAa4tLp8KzAKmAqcCN5pZza+quCv1bmHcdx98+tMD36f6ZyPFIlIsIsUivZQtjG8AP3L3\nLQDu/kZx+RnAPHff4u7LgR5gesJ6bCf1GAbsvoVx//1h6q2ISLNJmTAOBT5tZg+a2SIz+2hx+SSg\nfEi4t7hsSAxFC2PKFFixIoxXlHOHRx6B448f+D7VPxspFpFiESkW6VV1jVQzWwhMLF8EOHBlcd/j\n3P0EMzsO+CVwcDV/rxaGYgxj5EiYNAmWL4dDD43LV68OSUP3vRCRZlRVwnD3U3a1zsy+Dvy6uN0j\nZtZnZhMILYoDyzZtLy7bqe7ubjo6OgBoa2ujs7Nz2y+JUp/lQMovvQTQxahRg3t/peUpU+DXvy5w\nwglx/a23FjjwQAjDOwPbX3n/bIr6NlO5tKxR6lPP8pIlS7j00ksbpj71LM+ZM6fq74dmLRcKBebO\nnQuw7fsyhWS3aDWzvwcmuftsMzsUWOjuB5nZNOB24HhCV9RC4JCd3Yu1VrdoLXfTTXD++eECgfvs\nU9Ndb+eii+Cww+CSS+Kya66BlSvhuusGvr9CobDtg5J3ikWkWESKRZTqFq0pb9vzv4GbzewpYCNw\nLoC7LzWzXwBLgc3AhTXPCrsxFGMYAAcfDC++uP2ynp5wk6XB0H+ESLGIFItIsUgvWcJw983AObtY\ndzVwdaq/vTulMYxS4kilvR0efHD7ZT09A7vgoIhII8nlmd4jR6a/aVF7O/T2G5l5/nk45JDB7a+8\n/z7vFItIsYgUi/RymTBaEs6QKpk0KUytLXnvPXj99ZBIRESaUbJB71pIMej9u9/BV74Cf/5zTXe7\ng02bws2R3nsvXJW2pwdmzoQXXkj7d0VEUg16q4WRSEsLjB8Pr70Wyr29odUhItKscpkwUs+QKtl/\n/3CjJAjTaas5YU/9s5FiESkWkWKRXi4TxlC0MAAmTIhdXytXqoUhIs0tdwljxIiha2GUJ4ze3upa\nGJpjHikWkWIRKRbp5S5hTJwYzsAeCv1bGLqGlIg0s9wljIMPDvfUHgrlCeONN2C//Qa/L/XPRopF\npFhEikV6uUsYQ6k8Ybz1lu6yJyLNLXfnYQylW26BBQvgtttg8mT44x9DC0dEJCWdh9GExo+HN98M\nr9euVQtDRJqbEkZCpS6prVvh7bdh7NjB70v9s5FiESkWkWKRnhJGQmPGwF/+EpLFmDHhEiEiIs1K\nYxgJvfginHQSLFoEn/lMuGWriEhqGsNoQq2t8O67miElItmghJFQayusXx8GvMeNq25f6p+NFItI\nsYgUi/SUMBIaPTpc3vytt6ob8BYRaQQaw0hsr73ghhvCORh33FHv2ohIHmgMo0m1toYWxlBd8FBE\nJBUljMRKCaPaS6qrfzZSLCLFIlIs0lPCSKy1NZztrRaGiDQ7jWEkdtxx4fpR7e3wk5/UuzYikgca\nw2hSamGISFYoYSQ2erTGMGpNsYgUi0ixSE8JIzHNkhKRrFDCSKzUJVVtC0P3K44Ui0ixiBSL9JQw\nEmtthXXr1MIQkeanhJFYayu4awyjlhSLSLGIFIv0lDASa20Nz2phiEizU8JIbPTo8KwxjNpRLCLF\nIlIs0lPCSEwtDBHJimQJw8yOM7OHzezx4vOxZeuuMLMeM1tmZjNS1aERlBKGxjBqR7GIFItIsUhv\nRMJ9/ytwpbsvMLNTgWuAz5jZNGAWMBVoB+4xs0Oa/hogu6AWhohkRcouqVVA6bZBbUBv8fXpwDx3\n3+Luy4EeYHrCetRVrVoY6p+NFItIsYgUi/RStjAuBx4ws58ABny8uHwSsLhsu97iskwqDXqrhSEi\nza6qhGFmC4GJ5YsAB64ELgYudvc7zezzwM3AKQP9G93d3XR0dADQ1tZGZ2fntl8SpT7LRi4/+yxA\nFy0t1e2vvH+2kY6vHuXSskapTz3LS5Ys4dJLL22Y+tSzPGfOnKb7fqhVuVAoMHfuXIBt35cpJLu8\nuZm97e77lJXXunubmV0OuLv/uLj8LmC2uz+0k300/dDGsmUwbRo8/TQcccTg91MoFLZ9UPJOsYgU\ni0ixiFJd3jxlwvh/wLfd/T4zOxn4kbsfVxz0vh04ntAVtRDY6aB3FhLGK6/AQQfBc8/BIYfUuzYi\nkgepEkbKMYwLgP9uZi3ABuDvAdx9qZn9AlgKbAYubPqssBuaJSUiWZFslpS7P+rux7v70e7+MXd/\nvGzd1e4+xd2nuvuCVHVoBDoPo/YUi0ixiBSL9HSmd2KjRoGZWhgi0vx0T+8h8LWvwY03wsiR9a6J\niORB0w1610JWEoaIyFBKlTDUJdUk1D8bKRaRYhEpFukpYYiISEXUJSUikjHqkhIRkbpSwmgS6p+N\nFItIsYgUi/SUMEREpCIawxARyRiNYYiISF0pYTQJ9c9GikWkWESKRXpKGCIiUhGNYYiIZIzGMERE\npK6UMJqE+mcjxSJSLCLFIj0lDBERqYjGMEREMkZjGCIiUldKGE1C/bORYhEpFpFikZ4ShoiIVERj\nGCIiGaMxDBERqSsljCah/tlIsYgUi0ixSE8JQ0REKqIxDBGRjNEYhoiI1JUSRpNQ/2ykWESKRaRY\npKeEISIiFdEYhohIxmgMQ0RE6qqqhGFmnzezp82sz8yO6bfuCjPrMbNlZjajbPkxZvakmT1nZnOq\n+ft5ov7ZSLGIFItIsUiv2hbGU8B/Ae4rX2hmU4FZwFTgVOBGMys1j/4HcJ67HwocamafrbIOubBk\nyZJ6V6FhKBaRYhEpFulVlTDc/Vl37wH695WdAcxz9y3uvhzoAaab2f7A3u7+SHG7W4Azq6lDXqxd\nu7beVWgYikWkWESKRXqpxjAmAa+WlXuLyyYBK8qWryguExGRBjdiTxuY2UJgYvkiwIEfuPtvU1VM\ntrd8+fJ6V6FhKBaRYhEpFunVZFqtmS0CLnP3x4rlywF39x8Xy3cBs4GXgUXuPrW4/GzgRHf/xi72\nqzm1IiKDkGJa7R5bGANQXrn5wO1mdh2hy2kK8LC7u5mtM7PpwCPAucBPd7XDFAcsIiKDU+202jPN\n7FXgBOB3ZvZ/ANx9KfALYCnwB+DCsjPwLgJuAp4Detz9rmrqICIiQ6Ohz/QWEZHG0ZBnepvZTDP7\nU/Hkvu/Vuz6pmNlyM3vCzB43s4eLy8aZ2QIze9bM7jazsWXbZ+ZkSDO7yczWmNmTZctqduxm1mJm\n84rvWWxmBw7d0Q3MLmIx28xWmNljxcfMsnWZjIWZtZvZvWb2jJk9ZWaXFJfn7nOxk1hcXFxe38+F\nuzfUg5DEngcOAkYCS4DD612vRMf6IjCu37IfA/+1+Pp7wI+Kr6cBjxPGnTqKMSq1EB8Cjiu+/gPw\n2XofWwXH/kmgE3gyxbED3wBuLL7+IuG8oLof9wBiMRv49k62nZrVWAD7A53F12OAZ4HD8/i52E0s\n6vq5aMQWxnTC2MbL7r4ZmEc4ETCLjB1beWcA/158/e/EExtPJ0MnQ7r7/wXe6re4lsdevq9fASfX\n/CBqZBexgB1PiIXBnRTbFLFw99XuvqT4+i/AMqCdHH4udhGL0jlrdftcNGLC6H/SX5ZP7nNgoZk9\nYmbnF5dNdPc1ED40wH7F5Xk4GXK/Gh77tve4ex+w1szGp6t6Et80syVm9r/KumFyEQsz6yC0uh6k\ntv8nmjkWDxUX1e1z0YgJI08+4e7HAJ8DLjKzTxGSSLk8z0qo5bE32xTtG4GD3b0TWA38pIb7buhY\nmNkYwi/ebxV/Xaf8P9Fssajr56IRE0YvUD740l5cljnuvqr4/DpwJ6E7bo2ZTQQoNidfK27eC3yw\n7O2luOxqeTOq5bFvW2dmw4F93P3NdFWvLXd/3Yudy8D/JHw2IOOxMLMRhC/IW939N8XFufxc7CwW\n9f5cNGLCeASYYmYHmVkLcDbhRMBMMbPRxV8PmFkrMINw9d/5QHdxs68Apf8084GzizMbJhNPhlwN\nrDOz6WZmhJMhf0NzMHY84bO7+LraY59f3AfAF4B7kx1FbWwXi+IXY8lZwNPF11mPxc3AUne/vmxZ\nXj8XO8Si7p+Les8G2MUMgZmEWQE9wOX1rk+iY5xMmAH2OCFRXF5cPh64p3j8C4C2svdcQZj9sAyY\nUbb8o8V99ADX1/vYKjz+nwMrgY3AK8DfAeNqdezAKMLJoz2EfvCOeh/zAGNxC/Bk8TNyJ6EfP9Ox\nAD4B9JX9v3is+F1Qs/8TGYhFXT8XOnFPREQq0ohdUiIi0oCUMEREpCJKGCIiUhElDBERqYgShoiI\nVEQJQ0REKqKEISIiFVHCEBGRivx/fBdbd7eHj3MAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x106b9dc90>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "bb, aa = signal.ellip(6, 1, 40, 1300.0/24000, 'low')\n",
    "w,h = signal.freqz(bb, aa)\n",
    "w = w * 48000 / (2 *np.pi)\n",
    "pl.ylim((-100, 10))\n",
    "pl.plot(w, 20 * np.log10(abs(h)));\n",
    "bb,aa"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAD7CAYAAACIYvgKAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmYVNWd//H3VxYRUBsUsO0GGwQV3FqDZDExNRpxySgm\nk0Uni62JSR4dk0mMUUNG8mRM0EySiXGeJC5J2ixqlt8v7hutlk9+o4CJtqigtoEGaRExAoqKQHN+\nf5wqKdpearm3b91Tn9fz3Kfq3trOtw/Ut8753sWcc4iISG3bJekGiIhI8pQMREREyUBERJQMREQE\nJQMREUHJQEREgKFJN8DMtG+riEgZnHMW1XtVxcjAORfsMm/evMTboPgUXy3GF3JszkX/G7oqkkHI\nOjs7k25CrBRfuoUcX8ixxUHJQERElAzi1tLSknQTYqX40i3k+EKOLQ4Wx9xTSQ0wc0m3QUQkbcwM\nF1oBOWTZbDbpJsRK8aVbyPGFHFsclAxERETTRCIiaaRpIhERiZySQcxCn7dUfOkWcnwhxxaH2JOB\nmZ1oZk+b2bNmdlHcnyciIqWLtWZgZrsAzwLHAS8AjwCnO+eeLniOagYiIiVKW81gFtDhnFvpnNsK\n3ATMifkzRUSkRHEngwbg+YL11bltNSP0eUvFl24hxxdybHFI/BTW1aa7G155BV56Cdau9bcvvQQv\nvwwbNsDGjTuWV1+FzZv98tZbOxaAXXbxS3c37LorjBgBo0fDqFH+dvRo2HNPGDcO9t7b344bB+PH\nQ0MD7LsvDB+e7N+iFN3d8Nxz8Oyz0NGxY1mxwv8desZ+yCFw8skwc6Z/XESSFXcy6AImFaw35rbt\npKWlhaamJgDq6upobm4mk8kAO7J7uesPPJBl82aYMSPDSy9BW1uW9eth7NgMa9fCE0/49a1b/ePr\n1mUZPRoaGjJMmACQpa4Ompsz7L8/rFmTZepUOProDHvs4V8/fDgcc0yGESNg0aIsZv7x7dvhL3+B\n7u4sM2dm2LQJHnwwy5tvwtSpGTZsgIcfzvLUUzByZIZ166CjI5tLPBnq6mCPPbKMGweHH56hoQFe\ney3L3nvDySdnaGyExx7znxfV36uY9S1bYNy4DI8+CrfdBuedl2Xlygzjx8Pee2dpbIRjj81w2mn+\n7wlw8ME+/oceyvLGG7BpU4azzoKurixHHQVnn51h9mz/9+zr89esgUsvzfLww3DPPf7vEXe8+W2D\n+fcdzPWQ48tkMlXVnkrXs9ksra2tAG9/X0Yp7gLyEOAZfAF5DbAYOMM5t6zgOQMWkLu7yX2BwOuv\n+yV/f8MG/6u9r2XdOv8e+V/dEyb4277u7703DBsW11+keNu3+xHJ6tXQ1fXOJb+9u9uPJBoaoLFx\nx/3CbRMmwNAS0v7Wrf5vlx8VrVgBTz8Nzzzjb7u6YOpUOPJIOOIIf9vc7Ec6pershLvu8suDD8L+\n+8Nxx/nlAx/wI6p774Vrr4UHHoCPf9yPvoYNg+uuK/3zREIRdQE59iOQzexE4Ep8feIXzrnLezzu\nZs1ydHfz9rJtm//yz3/pv/UWjBy5Y6ohP90watTOUy09l/z2kSNjDbFfhb+64vDaa70nicJl7Vo/\n5VQ4TTNypE8427b5pbvbJ4H16/3011577UigkybB9Olw0EF+mTx5R8KMMr6tW2HxYrjvPr/87W++\nrRMnwjnnwBlnwO67+x8ABx7on3PIIZF8dJ/i7r+khRxfyLFB9Mkg9pqBc+5u4MD+nnPllTBkyM5L\n/kt/1CjYbTfNK/dl9913fEn3xTlf18iPpvK3Q4b4EUP+duhQGDMGxo5N5u89bBgcfbRfLr3Ut3HN\nGj8KKVRXB9/8Jlx0Edxxx+C3UyREOjeRpNKWLX60cu21cOyxSbdGZPCl7TgDkVgMHw7z58OFF/rp\nLhGpjJJBzPJ7A4Qqyfg+/nE/xXXTTfF9hvovvUKOLQ5KBpJaZvCDH/j6webNSbdGJN1UM5DUmzMH\njjkGLrgg6ZaIDJ7U7Vo6YAOUDKRCy5b5ZPDMM35PKJFaoAJyyoQ+b1kN8U2f7ndHvffe6N+7GuKL\nU8jxhRxbHJQMJAgHHOCPlBaR8miaSILws59BeztcfXXSLREZHJomEulFU5M/z5GIlEfJIGahz1tW\nS3xNTfFME1VLfHEJOb6QY4uDkoEEoakJVq3S0cgi5VLNQIKxzz7+TKcNNXUtPalVqhmI9GHyZNUN\nRMqlZBCz0Octqym+OOoG1RRfHEKOL+TY4qBkIMHQyECkfKoZSDCuucZfKU2Xw5RaoJqBSB/i2r1U\npBYoGcQs9HnLaoovjmmiaoovDiHHF3JscVAykGBMmgSrV0N3d9ItEUkf1QwkKI2N8NBDPjGIhEw1\nA5F+6BxFIuVRMohZ6POW1RZf1EXkaosvaiHHF3JscagoGZjZx8zsSTPrNrMjezx2iZl1mNkyM5td\nWTNFiqNjDUTKU1HNwMwOBLYDVwNfd849mts+HbgBOApoBNqAab0VB1QzkCj94hfwl79Aa2vSLRGJ\nV1XVDJxzzzjnOoCeDZoD3OSc2+ac6wQ6gFmVfJZIMTQyEClPXDWDBuD5gvWu3LaaE/q8ZbXFF3UB\nudrii1rI8YUcWxyGDvQEM1sATCjcBDhgrnPutiga0dLSQlNTEwB1dXU0NzeTyWSAHR2a1vX29vaq\nak/o8S1fnqWrC7ZuzTBsWHjxhd5/Wu97PZvN0pqb/8x/X0YpkuMMzOwB4IKCmsHFgHPOXZFbvxuY\n55xb1MtrVTOQSO23H2SzfspIJFRVVTPoobBRtwKnm9lwM5sMTAUWR/hZIn3SOYpESlfprqWnmdnz\nwHuA283sLgDn3FLgD8BS4E7g3Fr9+Z8f5oWqGuOLsohcjfFFKeT4Qo4tDgPWDPrjnLsZuLmPx+YD\n8yt5f5FyaGQgUjqdm0iCc/310NYGv/lN0i0RiU811wxEqoLOTyRSOiWDmIU+b1mN8UU5TVSN8UUp\n5PhCji0OSgYSnIYGWLcO3nor6ZaIpIdqBhKkKVPgnntg2rSkWyISD9UMRIqgcxSJlEbJIGahz1tW\na3xRFZGrNb6ohBxfyLHFQclAgjR5so41ECmFagYSpN/+Fu68E264IemWiMRDNQORIugoZJHSKBnE\nLPR5y2qNL6oCcrXGF5WQ4ws5tjgoGUiQ6uth/Xp4882kWyKSDqoZSLAOOABuvRUOOijplohETzUD\nkSLpHEUixVMyiFno85bVHF8UReRqji8KIccXcmxxUDKQYOkoZJHiqWYgwfrd7+D22+HGG5NuiUj0\nVDMQKdLEifD880m3QiQdlAxiFvq8ZTXHN2kSrFpV2XtUc3xRCDm+kGOLg5KBBKuhAV58EbZtS7ol\nItVPNQMJWkMDLFzop4xEQqKagUgJopgqEqkFSgYxC33estrjq7SIXO3xVSrk+EKOLQ4VJQMz+76Z\nLTOzdjP7P2a2R8Fjl5hZR+7x2ZU3VaR0GhmIFKeimoGZfQi43zm33cwuB5xz7hIzmwH8DjgKaATa\ngGm9FQdUM5A4XXklPPccXHVV0i0RiVZV1Qycc23Oue251YX4L36AU4GbnHPbnHOdQAcwq5LPEimH\nRgYixYmyZnA2cGfufgNQOFPbldtWc0Kft6z2+CpNBtUeX6VCji/k2OIwdKAnmNkCYELhJsABc51z\nt+WeMxfY6pwr68D/lpYWmpqaAKirq6O5uZlMJgPs6NC0rre3t1dVe2otvlWrsixfDhBmfKH3n9Z3\nrGezWVpbWwHe/r6MUsXHGZhZC3AOcKxz7q3ctovx9YMrcut3A/Occ4t6eb1qBhIb52DkSPjHP/yt\nSCiqqmZgZicCFwKn5hNBzq3A6WY23MwmA1OBxZV8lkg5zKCxUecoEhlIpTWDq4DRwAIze9TMfgrg\nnFsK/AFYiq8jnFurP//zw7xQpSG+SuoGaYivEiHHF3JscRiwZtAf59y0fh6bD8yv5P1FoqA9ikQG\npnMTSfAuvRR22QW+/e2kWyISnaqqGYikgUYGIgNTMohZ6POWaYhv0qTyC8hpiK8SIccXcmxxUDKQ\n4E2cqJGByEBUM5DgbdoE48fD66/7XU1FQqCagUiJRo+GESP8gWci0jslg5iFPm+ZlvjKLSKnJb5y\nhRxfyLHFQclAaoL2KBLpn2oGUhPOOw8OOgjOPz/plohEQzUDkTJoZCDSPyWDmIU+b5mW+Mq9FnJa\n4itXyPGFHFsclAykJmhkINI/1QykJqxaBe97H6xenXRLRKIRdc1AyUBqwrZt/uI2r78Ow4Yl3RqR\nyqmAnDKhz1umJb6hQ2HCBHjhhdJel5b4yhVyfCHHFgclA6kZ5RaRRWqBpomkZpx+Opx6Kvzrvybd\nEpHKaZpIpEzao0ikb0oGMQt93jJN8ZUzTZSm+MoRcnwhxxYHJQOpGRoZiPRNNQOpGY8+CmefDe3t\nSbdEpHKqGYiUSSMDkb4pGcQs9HnLNMW3116webO/8lmx0hRfOUKOL+TY4lBRMjCz75jZ42bWbmZt\nZtZY8NglZtZhZsvMbHblTRWpjJmONRDpS0U1AzMb7ZzblLt/PnCYc+4cM5sB/A44CmgE2oBpvRUH\nVDOQwXT88fD1r8MJJyTdEpHKVFXNIJ8IckYB+avMngrc5Jzb5pzrBDqAWZV8lkgUJk5U3UCkNxXX\nDMzsMjNbBbQA83ObG4DCwXhXblvNCX3eMm3xTZpU2jRR2uIrVcjxhRxbHIYO9AQzWwBMKNwEOGCu\nc+4259y3gG+Z2UXAj4GzSm1ES0sLTU1NANTV1dHc3EwmkwF2dGha19tz+zFWS3tqPb433siyZAlA\nmPGF3n+1vJ7NZmltbQV4+/sySpEdZ2BmE4E7nXOHmtnFgHPOXZF77G5gnnNuUS+vU81ABs2CBXD5\n5XDffUm3RKQyVVUzMLOpBaunAfnDeW4FTjez4WY2GZgKLK7ks0SiMH48rFuXdCtEqk+lNYPLzWyJ\nmT2GH3dfAOCcWwr8AVgK3AmcW6s///PDvFClLb6xY+GVV4p/ftriK1XI8YUcWxwGrBn0xzn3sX4e\nm8+OgrJIVSg1GYjUCp2bSGqKczBiBGzYALvtlnRrRMpXVTUDkbQx86elWL8+6ZaIVBclg5iFPm+Z\nxvjGjoV//GPg50E64ytFyPGFHFsclAyk5qhuIPJOqhlIzTntNDjzTPjIR5JuiUj5VDMQqZBGBiLv\npGQQs9DnLdMY3157qWaQF3J8IccWByUDqTkaGYi8k2oGUnOuvhr+9je45pqkWyJSPtUMRCqkkYHI\nOykZxCz0ecs0xrfXXsUngzTGV4qQ4ws5tjgoGUjNKeWgM5FaoZqB1JxVq+Doo0u74plItVHNQKRC\nqhmIvJOSQcxCn7dMY3yjRsG2bbB588DPTWN8pQg5vpBji4OSgdQcM40ORHpSzUBq0sEHw+9/D4cc\nknRLRMqjmoFIBDQyENmZkkHMQp+3TGt8xR5rkNb4ihVyfCHHFgclA6lJGhmI7Ew1A6lJX/86TJgA\nF16YdEtEyqOagUgENDIQ2ZmSQcxCn7dMa3zFJoO0xleskOMLObY4RJIMzOwCM9tuZmMLtl1iZh1m\ntszMZkfxOSJRKeVkdSK1oOKagZk1AtcBBwLvcs69YmbTgRuAo4BGoA2Y1ltxQDUDScJ998F3vwv3\n3590S0TKU401g/8Gepbh5gA3Oee2Oec6gQ5gVgSfJRIJ1QxEdja0kheb2anA8865J8x2SlANwMMF\n6125bTUnm82SyWSSbkZs0hpfKTWDNMZXrGLie+MNePFFv2za5M/pVLg4B8OGwfDh/nbYMNhtN6ir\n23kZMcKfCmSwhN53URswGZjZAmBC4SbAAd8CvgkcX2kjWlpaaGpqAqCuro7m5ua3OzFfBErrent7\ne1W1R/H59ZkzM7zySv/Pf+45uPvudrZsgdmzq6v9UfffzJkZnnwS/vjHLMuXw8aNGbq6YPXqLFu3\nwr77ZthnH9i6Ncvw4dDYmGHECHjlFf9+48Zl2LoVurqybNsGo0Zl2LgRXnghy6ZN8OabGZyDuros\nY8fCgQdmqK+HzZuzTJgAJ52UoakJOjqyDBlSPX+falrPZrO0trYCvP19GaWyawZmdgi+FvAGPkE0\n4kcAs4CzAZxzl+eeezcwzzm3qJf3Uc1ABp1zsOuu8Npr/rY3RxzhL4Kzdi3ssw9Mm+aXc8+FQw8d\n3PZG7bXXoK0N7rgDsll44QWYMQMOO8wvBx8MjY1QXw977hnNL/o33vB/yzVrdl5WroTOTlixAl5+\n2X/u/vvDgQfuvDQ2wi7a//FtUdcMIjvozMxWAEc659ab2Qzgd8C78dNDC1ABWarMPvvAY4/5L7ze\n1NX5L6jdd/dfWB0d/ovznnvgr3+FIUMGtbkV6+qCP/7RJ4CFC+E974EPfxiOP95/2Q6taNI4Gm+9\n5f/Wzz0Hzzzjl2ef9bevvuqT8OGH71gOPRRGj0661cmoxgJynsOPEHDOLQX+ACwF7gTOrdVv/Pww\nL1Rpjq+/usGGDX700N6eZehQ/0v1xBNh/nwYORJ+/evBbWslHnsMPv1p/8X55JNw3nl+JLBgATQ3\nZzn44OpIBOBHaQccACefDF/9Kvz8536Pr64uf4W6+fNh+nR45BE4/3wYP96P1j72MfjP/4Rbb/XP\ncy7d/zaTENk/AefclB7r84H5Ub2/SNT6O9agsxOamt45PWIGP/oRfPSj8PGPV++v0u3b4e674Yc/\n9L+qv/xl+J//8aOdtBozBj74Qb/kbdvmRw6PP+6Xn/3MJz+AqVPhn/8Z3vtemDnTX9RI+qZzE0nN\nmjMHzj7b3/Z0883wq1/BLbf0/tpPfcqPFr7znXjbWI6nnoIzz/RflBdcAJ/8pN/Tp1Y4569v/fDD\nfjps4UJYssRPhb3nPX5573t9shjMvZuiFvU0UZUMDkUG39ixvkDcm85O2G+/vl87f74vMJ9zDkyc\nGEvzStbdDT/+MVx+uW/f5z6X7i+7cpnBpEl++eQn/ba33vIjhoULfc3kP/4DXn8d3vc++Kd/gkzG\n1yBquUBdw6EPjtDnLdMcX381g/w0UV/xTZoEX/oSzJ0bV+tKs2KF/1K75RZYvBg+//niEkGa+28g\nhbHtuqsfEfz7v/sr3K1c6UcLn/qUn2Y64wwYNw4+8hG48kr/2PbtybU9CUoGUrOKSQb9ufhiX4T9\n61+jbllpfvUrmDXLT3c98ABMnpxse9Ji3339yOFnP4Onn4YnnoBPfMJPs33sY744/S//4mstTz3l\np59CppqB1Kyf/xza2/1tT83N8MtfwpFH9v8e114Lv/kNPPhgMlMyP/kJXHWVr3EcfPDgf37IVq/2\nuxJns/5cVlu3+j3KTjoJjjsu+WJ8Ne9aKpIq/dUMVq4ceGQAvgC9YQP8+c+RNq0ov/41/OAH/uAx\nJYLoNTb6XXKvuw6WL/e7uB56qF+fOBE+8AH43vfg0UfDmFJSMohZyHOykO74+pom2rDBF2PHjBk4\nviFD4LLL/JfyYLr5ZrjoIn8AXH+F7oGkuf8GEmVsZv74h698Be66C156ydeLXnrJ1xv23dfvwXXj\njek9AaKSgdSsvpJBflRQ7LTP7Nn+YK6XX460eX267z74whfg9tv9AVgy+HbbzU8Z/fjH/jiOhx6C\nd78bbrjB/9s57jg/fff880m3tHiqGUjNWrkSjjnG3xa65Rb4xS/80azFmjPHH4T26U9H28aeFi2C\nU06BP/3Jt12qz+uv+x0L/vxnn7AnT4bTTvN7Ks2YEV1tSTUDkYj0NTIoZk+inj78Ybjzziha1bdV\nq3zSaW1VIqhmo0b5L//rr/cn5vuv//LTSSed5KeaLrzQjySqrc6gZBCzkOdkId3xjR7tz8e/ZcvO\n2wuTQbHxnXyyn7/fti3KFu7skkvgi1/0nxWVNPffQKohtqFD/fEfP/mJH4H+/vd+iumLX4SGBn+a\nkIcfro7dVpUMpGaZ9T46KGdk0Njo9zBZuDCq1u3skUf8Lo4X9rymoKSGmd9V+Tvf8cc0PPgg7L23\n3yNtyhSf7JcsSS4xqGYgNW3GDD//PmPGjm1HHOFrBgMdY9DT3Ll+6D8/4tMzOudPznbmmf4UExIW\n5/xJ9m68EW66yU8znXGGX6ZO7ft1qhmIRCiqkQH4usEdd0TRqp3dcgusXw8tLdG/tyTPzB/keMUV\n/rQi113nawzvfz8cdZQ/S+6aNfG3Q8kgZtUwbxmntMfX88CzwmMMoLT43v1uf52AVauia9+WLfCN\nb/jjGOK4mE7a+68/aYxtl138yfOuusofAT1/vt9tecYMOPVU/8Ng69aYPjuetxVJh54jg5Ur/UFc\n5ez+N2SI3/c8yr2Krr7azyefcEJ07ynpMHQofOhD/rQozz/vd039wQ98beob34j+81QzkJr2ta/5\nvTouuMCv33KLH6bfdlt573fjjf7Ao3JfX2jDBn8O/ra29F9zWaLzzDM+QXz/+6oZiESm59XOyq0X\n5J1wgt9L5M03K22ZP+/NqacqEcjODjzQ1xeipmQQszTOW5Yi7fH1nCbqmQxKjW/sWH+RlEr/LJ2d\nfo+muK+klvb+60/IscVByUBqWs8CcqUjA4hmr6LWVvjsZ6G+vrL3ESmWkkHMMplM0k2IVdrjG2hk\nUE58+WRQSSmsrS3aI437kvb+60/IscVByUBqWtQ1A4BDDvG7py5bVt7rX33VH4T0/vdX1g6RUigZ\nxCz0ecu0x1c4Mti40Z9baOzYHY+XE59ZZVNFDz7oj1nYbbfyXl+KtPdff0KOLQ4VJQMzm2dmq83s\n0dxyYsFjl5hZh5ktM7PZlTdVJHqFNYNSr2PQn0qSwYIFfv9ykcFU0XEGZjYPeM0596Me26cDNwBH\nAY1AGzCttwMKdJyBJMk5GD7cn4P+7rv9NY2jOEZg40Z/9atNm0pPLjNm+Osqv+tdlbdDwlWN5ybq\nrTFzgJucc9ucc51ABzArgs8SiVT+zKXr10dTL8jbc0//3q+9Vtrrurr8eWmam6Nph0ixokgG/2Zm\n7WZ2nZntmdvWABRe8K0rt63mhD5vGUJ8+bpBb8mgkvjq60s/wVhbGxx7bDznIepNCP3Xl5Bji8PQ\ngZ5gZguACYWbAAfMBX4KfMc558zsMuCHwOdLbURLSwtNuf+FdXV1NDc3v71bWL5D07re3t5eVe1R\nfO9cHzIEXnklQ2cn7Llnlmw2mvjq6+Guu7KsWVP863/72yyHHQYwOPGH0H+1sp7NZmltbQV4+/sy\nSpGdm8jM9gNuc84dZmYXA845d0XusbuBec65Rb28TjUDSdQpp8A558C3vw3XXAMzZ0bzvp/8pL/8\n4RlnFPd853yd4aGH/HVzRfpTVTUDM9unYPWjwJO5+7cCp5vZcDObDEwFFlfyWSJx6W+aqBL19f6U\n1sV66ikYOVKJQJJRac3g+2a2xMzagQ8CXwVwzi0F/gAsBe4Ezq3Vn//5YV6oQohvr738RUW2bPH3\nC1USX6k1gyR2KQ2h//oScmxxGLBm0B/n3Gf7eWw+EPEFAEWiN3YsLF4c3TEGefX1/pq2xWprg7PO\niu7zRUqh6xlIzfvpT+G73/XXPr799ujed8ECf6Wq++8f+LlbtviLo69Y8c7RiUhvqqpmIBKCsWP9\n3H7UO2jsu2/x00QLF/rz1CsRSFKUDGIW+rxlCPHlv4B7SwaDVTNI6hQUIfRfX0KOLQ5KBlLz8iem\ni3pkMGYMbN5c3FXP2trg+OOj/XyRUqhmIDVvxQp/0flHHonuGIO8piZfM5gype/nbNwIjY2wbh2M\nGBHt50u4VDMQiVhcIwMo7liDBx6A971PiUCSpWQQs9DnLUOIb489YO7c3ou3lcZXTN3g/vvhuOMq\n+piyhdB/fQk5tjgoGUjNM4PLLov2GIO8YpLBsmVw+OHRf7ZIKVQzEInRZZf5ayXM7+fwyylT4J57\nYNq0wWuXpJ9qBiIpMtCxBlu3+msY7Lff4LVJpDdKBjELfd5S8fVvoGmiVav8c4YPr+hjyhZy/4Uc\nWxyUDERiNFAyWL4c9t9/8Noj0hfVDERitHYtHHKIP4agNz//Ofztb/7ayyKlUM1AJEXGjYMNG/yJ\n6HqjkYFUCyWDmIU+b6n4+rfLLjB+PLz4Yu+P//3vySaDkPsv5NjioGQgErP+6gZ//3v/p6oQGSyq\nGYjE7JRT4HOf89dDLuQc7Lmn36Oori6Ztkl6qWYgkjJ9HWvw8sswdKgSgVQHJYOYhT5vqfgG1tc0\nUTUUj0Puv5Bji4OSgUjM+koGSRePRQqpZiASs1tvhWuueef1lS+7DN54A773vWTaJemmmoFIyvR1\nTQONDKSaVJwMzOx8M1tmZk+Y2eUF2y8xs47cY7Mr/Zy0Cn3eUvENrL9poqR3Kw25/0KOLQ5DK3mx\nmWWAU4BDnXPbzGzv3PbpwCeA6UAj0GZm0zQfJLVowgS/51B3NwwZsmN7NRSQRfIqqhmY2e+Bq51z\n9/fYfjHgnHNX5NbvAr7tnFvUy3soR0jwxo+Hxx/3owSAN9+EMWP8tQ4KE4RIsaqtZnAAcIyZLTSz\nB8zsXbntDcDzBc/rym0TqUk9p4o6O/01DJQIpFoMmAzMbIGZLSlYnsjdnoqfZhrjnHsP8A3gj3E3\nOG1Cn7dUfMXpeeBZtRSPQ+6/kGOLw4A1A+fc8X09ZmZfAv5v7nmPmFm3me2FHwlMKnhqY25br1pa\nWmhqagKgrq6O5uZmMpkMsKND07re3t5eVe1RfMnEV1+fYc2aHevLl2eYMiWc+LQe/3o2m6W1tRXg\n7e/LKFVaM/gC0OCcm2dmBwALnHP7mdkM4HfAu/HTQwuAXgvIqhlILZg7F3bdFS691K9/5SvQ1ARf\n/WqizZIUq7aawa+AKWb2BHAD8FkA59xS4A/AUuBO4Fx940st61kzqIbdSkUKVZQMnHNbnXOfcc4d\n6pyb6Zx7sOCx+c65qc656c65eytvajrlh3mhUnzF6ZkMqmW30pD7L+TY4qAjkEUGQWEy2L4dVqzQ\nyECqi85NJDIIVqyAD37QX7ugqwtmzuz7gjcixai2moGIFKG+3l/60rnq2a1UpJCSQcxCn7dUfMUZ\nMQJGj4ZOFzOVAAAEjUlEQVR//MPXC6pliijk/gs5tjgoGYgMknzdQCMDqUaqGYgMkg99CC68EK6/\nHk46CT7zmaRbJGmmmoFISuVHBtWyW6lIISWDmIU+b6n4ileN00Qh91/IscVBySBm+XO/hErxFa++\nHp591p++evz4yN62IiH3X8ixxUHJIGYbNmxIugmxUnzFq6+H//1fvyeRRTbTW5mQ+y/k2OKgZCAy\nSOrroaOjenYrFSmkZBCzzs7OpJsQK8VXvH339bfVUi+AsPsv5NjiUBW7libaABGRlIpy19LEk4GI\niCRP00QiIqJkICIiCScDMzvRzJ42s2fN7KIk21IuM+s0s8fN7DEzW5zbNsbM7jWzZ8zsHjPbs+D5\nl5hZh5ktM7PZybW8d2b2CzNba2ZLCraVHI+ZHWlmS3J9++PBjqMvfcQ3z8xWm9mjueXEgsfSFl+j\nmd1vZk+Z2RNm9uXc9tT3YS+xnZ/bHkT/mdmuZrYo913ylJl9L7d9cPrOOZfIgk9EzwH7AcOAduCg\npNpTQRzLgTE9tl0BfCN3/yLg8tz9GcBjwFCgKRe/JR1Dj7a/H2gGllQSD7AIOCp3/07ghKRj6ye+\necDXennu9BTGtw/QnLs/GngGOCiEPuwntpD6b2TudgiwEDh6sPouyZHBLKDDObfSObcVuAmYk2B7\nymW8c4Q1B7g+d/964LTc/VOBm5xz25xznUAH/u9QNZxz/w9Y32NzSfGY2T7A7s65R3LP+3XBaxLV\nR3zg+7GnOaQvvhedc+25+5uAZUAjAfRhH7E15B4Opf/eyN3dFf+9sp5B6rskk0ED8HzB+mp2dGya\nOGCBmT1iZp/PbZvgnFsL/h8wkD/5QM+Yu0hHzONLjKcB3595aejbfzOzdjO7rmAYnur4zKwJPwpa\nSOn/Jqs6xoLYFuU2BdF/ZraLmT0GvAhknXNLGaS+UwG5ckc7544ETgbOM7MP4BNEodD23w0tnp8C\nU5xzzfj/hD9MuD0VM7PRwJ+Ar+R+RQfzb7KX2ILpP+fcdufcEfjR3AfMLMMg9V2SyaALmFSw3pjb\nlirOuTW523XAzfhpn7VmNgEgN2R7Kff0LmBiwcvTEnOp8aQqTufcOpebXAWuZcfUXSrjM7Oh+C/L\n3zjnbsltDqIPe4sttP4DcM69ip/rn8kg9V2SyeARYKqZ7Wdmw4HTgVsTbE/JzGxk7lcKZjYKmA08\ngY+jJfe0M4H8f8hbgdPNbLiZTQamAosHtdHFMXaegy0pntxQdqOZzTIzAz5b8JpqsFN8uf9geR8F\nnszdT2t8vwSWOueuLNgWSh++I7ZQ+s/M9s5PcZnZbsDx+ALx4PRdwpXzE/F7BHQAFyfZljLbPxm/\nF9Rj+CRwcW77WKAtF9u9QF3Bay7BV/2XAbOTjqGXmG4AXgDeAlYBZwFjSo0HeFfub9IBXJl0XAPE\n92tgSa4vb8bP0aY1vqOB7oJ/l4/m/p+V/G+y2mLsJ7Yg+g84NBfTY8DjwNdz2wel73Q6ChERUQFZ\nRESUDEREBCUDERFByUBERFAyEBERlAxERAQlAxERQclARESA/w8keIKrwprxAQAAAABJRU5ErkJg\ngg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x1069ffd90>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "pl.xlim(0, 3000)\n",
    "pl.ylim((-65, 5))\n",
    "pl.plot(w, 20 * np.log10(abs(h)));"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([ 0.00107669, -0.00578907,  0.01354389, -0.0176606 ,  0.01354389,\n",
       "        -0.00578907,  0.00107669]),\n",
       " array([  1.        ,  -5.79629102,  14.04709981, -18.2174942 ,\n",
       "         13.33384205,  -5.22208148,   0.85492756]))"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAD7CAYAAAB5aaOHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XuUVOWZ7/Hv0zdic2sEJQrRFsEjELW9gdEk1mhEzMRo\nEkfJchJbc04uGhONM1HUjFkm42WSc4zJijGexMXEy5CLJpqsRMVIxThewAgHEJAGg3JR8QIoNDQN\nvOePt4q3KLqbpqt27121f5+1alW9e+/a9dZDUU+/z7v3LnPOISIisjc1cXdAREQqgxKGiIj0ihKG\niIj0ihKGiIj0ihKGiIj0ihKGiIj0Sl3cHeiJmemYXxGRPnDOWbn3mfgRhnNON+e44YYbYu9DUm6K\nhWKhWPR8i0riE4Z4K1eujLsLiaFYBIpFoFhETwlDRER6RQmjQrS2tsbdhcRQLALFIlAsomdR1rtK\nZWYuyf0TEUkiM8OlcdJbvGw2G3cXEkOxCBSLQLGInhKGiIj0ikpSIiJVRiUpERGJlRJGhVB9NlAs\nAsUiUCyip4QhIiK9ojkMEZEqozkMERGJlRJGhVB9NlAsAsUiUCyiF1vCMLOpZrbUzJaZ2dVx9UNE\nRHonljkMM6sBlgGnA2uBucA059zSou00hyEiso+qbQ5jEtDmnHvFOdcJzATOiakvIiLSC3EljFHA\nqoL26twy6Ybqs4FiESgWgWIRvYqZ9F60CK69Fr7xDXjmma63+dWv4Dvfgc2bffuuu6C9vf/6KCJS\nzeL6Te81wCEF7dG5ZXv4zGdaWbGimbY2OPnkJg49tIULLsjwoQ/B+ednGT4cTjklw5VXwkMPZRkz\nBp58MsOPfgRf+lKWzZvhyiszQPgLJJOpvHYmk0lUf9ROTjsvKf2Jq51flpT+9Gc7m80yY8YMAJqb\nm4lKXJPetcBL+Env14A5wGedc0uKtnPDhzsuvRSuugqGDvXLt2yB734X7rwTpk6F+fNhzBi4914Y\nNAjOPBPmzYN33oG774aLL+7nNygiEqOoJr1jO9PbzKYCt+PLYj93zt3SxTZu40bHkCFd7+OVV+CJ\nJ+Cww+DUU8Fy4dmyBR5/3CeNrVvhppsiexv9pvAvp7RTLALFIlAsgqgSRlwlKZxzjwD/Y2/bdZcs\nAA49tOvRw377wdlnw6ZN8OCDJXRSRER2qeprST3zDFxxBTz3XBk7JSKScNV2Hka/GDIE3nsv7l6I\niFSHqk4YgwZVT8IoPiImzRSLQLEIFIvoVXXCGDzYz2OIiEjpqnoOY9s2aGyEzs5wBJWISLXTHEYf\nNDRAba0/tFZEREpT1QkDqqcspfpsoFgEikWgWESv6hNGNU18i4jEqarnMACOOgruuw+OPrpMnRIR\nSTjNYfSRRhgiIuVR9QlDcxjVR7EIFItAsYheKhKGRhgiIqWr+jmMiy6CTEaXOBeR9NAcRh9VS0lK\nRCRuVZ8wqmXSW/XZQLEIFItAsYhe1ScMjTBERMqj6ucwfvhDaGuDH/2oTJ0SEUk4zWH0UbWUpERE\n4lb1CaNaSlKqzwaKRaBYBIpF9Ko+YWiEISJSHlU/h/HUU3D11fDf/12mTomIJJzmMPpIZ3qLiJRH\n1SeMailJqT4bKBaBYhEoFtGr+oRRLZPeIiJxq/o5jPZ22H9//UyriKSH5jD6aL/9YPt26OyMuyci\nIpWt6hOGmZ/HqPSylOqzgWIRKBaBYhG9qk8YUD0T3yIicar6OQyA8ePhgQdgwoQydEpEJOE0h1EC\njTBEREqXioRRDYfWqj4bKBaBYhEoFtFLTcLQCENEpDSRzWGY2X8AZwMdwArgYufcu7l104FLgO3A\n151zj3Wzj7LMYVx4IUydCp/7XMm7EhFJvEqcw3gMmOicawHagOkAZjYBOB8YD5wF3GFmZX9jhaqh\nJCUiErfIEoZz7nHn3M5c81lgdO7xJ4GZzrntzrmV+GQyKap+QHVMeqs+GygWgWIRKBbR6685jEuA\nP+YejwJWFaxbk1sWGY0wRERKV1fKk81sFjCycBHggOucc7/PbXMd0Omc+6++vEZrayvNzc0ANDU1\n0dLSQiaTAcJfFHtrDx6cYdWq3m+fxHYmk0lUf9ROTjsvKf2Jq51flpT+9Gc7m80yY8YMgF3fl1GI\n9MQ9M2sF/hdwmnOuI7fsGsA5527NtR8BbnDOPdfF88sy6X3XXTBnDvzsZyXvSkQk8Spu0tvMpgL/\nCnwynyxyHgammVmDmR0GjAXmRNUPqI7Daov/mkwzxSJQLALFInollaT24kdAAzArdxDUs865S51z\ni83sV8BioBO4tCzDiB4MGVL5CUNEJG6puJbUX/8K06f73/cWEal2FVeSSpIhQ+Ddd+PuhYhIZUtF\nwhg6tPIThuqzgWIRKBaBYhG9VCQMjTBEREqXijmMzk7/U62dnf4X+EREqpnmMEpQXw8NDdDeHndP\nREQqVyoSBlR+WUr12UCxCBSLQLGInhKGiIj0SirmMABOPBHuuMPfi4hUM81hlEgjDBGR0ihhVAjV\nZwPFIlAsAsUieqlKGBs3xt0LEZHKlZo5jMsvh3Hj4GtfK8vuREQSS3MYJar0kpSISNyUMCqE6rOB\nYhEoFoFiEb3UJIxquAChiEicUjOHcf/98Ic/+HsRkWqmOYwSVXpJSkQkbqlKGJV8WK3qs4FiESgW\ngWIRvVQlDI0wRET6LjVzGC+/DKefDn//e1l2JyKSWJrDKJGOkhIRKU1qEsbgwT5hJHhA1SPVZwPF\nIlAsAsUieqlJGA0NUFcHW7fG3RMRkcqUmjkMgJEjYcECfy8iUq00h1EGlX5orYhInFKXMCp14lv1\n2UCxCBSLQLGIXqoSho6UEhHpu1TNYZx7LrS2+nsRkWqlOYwyqOSSlIhI3JQwKoTqs4FiESgWgWIR\nPSUMERHplVTNYdxyC6xfD7feWrZdiogkTsXOYZjZVWa208z2L1g23czazGyJmU2Jug95OkpKRKTv\nIk0YZjYaOAN4pWDZeOB8YDxwFnCHmZU9E3alkktSqs8GikWgWASKRfSiHmHcBvxr0bJzgJnOue3O\nuZVAGzAp4n4AlZ0wRETiFlnCMLNPAquccwuLVo0CVhW01+SWRa6SE0Ymk4m7C4mhWASKRaBYRK+u\nlCeb2Syg8FJ+BjjgeuBafDmqJK2trTQ3NwPQ1NRES0vLrg9Gfgja2/ayZVnWrAHo2/PVVltttZPY\nzmazzJgxA2DX92UUIjlKysw+CDwOtOOTyGj8SGIScAmAc+6W3LaPADc4557rYj9lPUpqxQqYMsXf\nV5psNrvrg5J2ikWgWASKRVBRR0k55xY5597vnBvjnDsMWA0c65xbBzwMXGBmDWZ2GDAWmBNFP4rp\narUiIn3XL+dhmNnLwAnOuXdy7enAF4BO4OvOuce6eV5ZRxgdHf6X9zo6oH+OyxIR6X9RjTBSdeIe\nwPveBxs2+HsRkWpUUSWpJKvUI6XyE1yiWBRSLALFInpKGCIi0iupK0kddxz87Gf+XkSkGqkkVSYa\nYYiI9E3qEsbQoZV5aK3qs4FiESgWgWIRvdQlDI0wRET6JnVzGJddBhMm+HsRkWqkOYwy0QhDRKRv\nlDAqhOqzgWIRKBaBYhE9JQwREemV1M1h3HsvPPoo3HNPWXcrIpIYmsMoE40wRET6JpUJQ+dhVDbF\nIlAsAsUieqlLGIMHw3vvxd0LEZHKk7o5jMWL4bzz/L2ISDXSHEaZ7LcftLfH3QsRkcqTyoSxZUvc\nvdh3qs8GikWgWASKRfSUMEREpFdSN4exbRsMHAidnWXdrYhIYmgOo0zq62HnTti+Pe6eiIhUltQl\nDLPKLEupPhsoFoFiESgW0UtdwoDKTBgiInFL3RwGwKGHwpNP+nsRkWqjOYwy0ghDRGTfKWFUCNVn\nA8UiUCwCxSJ6qU0YOttbRGTfpHIO4/TTYfp0+NjHyr5rEZHYaQ6jjCqxJCUiEjcljAqh+mygWASK\nRaBYRE8JQ0REeiWVcxhf+hIceyx8+ctl37WISOwqcg7DzC43syVmttDMbilYPt3M2nLrpkTZh640\nNPiLEIqISO9FljDMLAOcDRzlnDsK+H5u+XjgfGA8cBZwh5mVPRP2pLYWduzoz1csneqzgWIRKBaB\nYhG9KEcYXwFucc5tB3DOvZVbfg4w0zm33Tm3EmgDJkXYjz3U1vor1oqISO9FmTCOAD5qZs+a2Wwz\nOz63fBSwqmC7Nbll/aampvJGGJlMJu4uJIZiESgWgWIRvbpSnmxms4CRhYsAB1yf2/cw59xJZnYi\n8GtgTCmvVy75ktSbb8LnPw/33gvDh+/9ec7BrFnw85/DkiX+Uuljx8KUKfDZz8KQIdH3XSQKzvnf\niNm+3f+4WP6+8HFPy7pat2OHH8nv2LHn457W9Xa7ve3Dua5vca4rjnlP7bo6WLYs2n/3fVVSwnDO\nndHdOjP7MvBgbru5ZrbDzIbjRxSHFGw6OresS62trTQ3NwPQ1NRES0vLrr8k8jXLfW3X1mbYuROu\nuy5LNgs33pjh9tt7fn5HB5x7bpYXX4Trr88wfTo8/3yWv/8dHn/ct884I8vnPgef+ERp/euqXVif\nLcf+KrmdX5aU/vS1/cQTWTo64IQTMrS3+/Vbt8KECRm2bIE5c/z6ww/3n79Fi7J0dsLo0b7d1pZl\n7dr5TJhwBR0d8Oqrfv2QIX79unW+PWCAb2/Y4Ns1NRk6O2Hr1izbt8POnRl27ICamiy1tX57/0Nj\nWerqYNCgDHV1sG2bXz9smG+3t/v1I0b47Tds8O2DD85QWwtvvZWlpgZGjfLt11/37UMO8e21a7OY\nwWGH+farr/r1Y8dmqKmBlSv96x1xhG+//LJff+SRfvu2Nt+eONG3H3nkB4wZ08LRR/vtFy3y+z/6\n6AxmsHCh3/6YY3x7wQK/vqUltAGOO86358/37eOP9/t74QW//Qkn+PV/+5tvT5rk288/v3t77lz/\nepMn+/acOX795Mn+3/+55/z+TzrJt5991q8vbGezvf9+mDFjBsCu78soRHZYrZl9ERjlnLvBzI4A\nZjnnDjWzCcB9wGR8KWoWMK6r42ejOqz2W9/yv7z31FNw6qlwzz2wdGn32zsH06ZBRwfcd5//iddi\nq1fDd74DDz0E3/se/PM/+xFIuWSzWQ25c+KIxbZt8N57/vbuu72737TJX7NsyxZ/X3zbts2fE9TY\n2P3tfe+DAQP2vDU0+PtVq7JMnJjpcpvibQvbdXX+/0Dhff8eelJ++j8SRHVYbZQJox64G2gBOoCr\nnHN/ya2bDnwB6AS+7px7rJt9RJIwvv1tnwR++1uYMQPOOAMWLIBR3cyk3HsvfP/78Oyz/j9wT55/\nHi66CI45Bu68U2WqpNm0Cd54A95+G955x9/Wrw+Pi9vr1/vbjh0weLD/9yy872nZwIH+1l0yGDCg\n8r+kJZmiShgllaR64pzrBD7XzbqbgZujeu29qa31f9298QYcdBCceCK88ELXCaOjA669Fn75y70n\nC4ATToC5c+HKK2HyZPjDH+Dww8v/HiTo7IS1a/0ob+1a/+/6xhuwbl14nL85ByNHwogRsP/+u98+\n8AGf6PffH4YNC8ubmvy/vb7cJe0iSxhJVlPjE8Y778ABB8DEifDii3D22Xtu+8tfwpFHwoc+1Pv9\nNzbCT38KP/kJnHIKPPggnHxyaX1O83B782ZYsQKWL4dXXoGnn84CGVatglWr/MELI0f6L/yDD4YD\nD/TtY4/19/nbgQfCoEHV9cWf5s9FMcUieqlMGLW1/q/PpiZfu504Ebo75+fuu+Hyy/v2Ol/5CjQ3\nwznnwG9+4+dLpGvbt/uEsGgRvPSSf5xPEuvXw5gxfqTW3OyPaMtkfILIJ4m6VH6SRfpXKq8l9b3v\nwZ//7EsYixb5uYmvftXPPxRatw7GjfP3Awb0/fWeeAIuuABmzvS/xZF269b5WC9c6G+LFvnDBw8+\nGD74QT+iGzs23A4+2I8KRaR3Km4OI8lqanyte2TuDJIxY2Dlyj23+9Of/I8slZIsAE47zZelPv1p\nePjhfStvVbr2dvjb32DOnHDbsMHP9RxzjI/vFVfA+PFdH30mIsmRyr/bamt9mWPoUN8+4ABfJ9+8\nefftnnjCH0FVDh/5CPziF/CpT/XtZJxKuU5ORwc8+aQ/Eu2jH/XzBv/yL37u4eyz4dFH/RFKs2b5\nI89aW33y2JdkUSmx6A+KRaBYRC+VI4z8UVL5ureZr4WvWuXLIXnPPgtXXVW+1z3rLLjpJn//9NNh\nhFPpXnsNfv97+N3v4K9/9TE87TS47jr48Ic1chCpFqlMGDU1/i/h2tqw7JBD4NVXQ8J4+23/RThx\nYnlf+5JLfPnrvPP8CKa+vnfPS9rRH6tXw/33+3NZli71SfCii/yJjcOGRfvaSYtFnBSLQLGIXmpL\nUoUjDPAJ45VXQnvBAjj66N2TSrl8+9u+HPbNb5Z/31Fqb/dJYsoUP/+wfDnceKM/v+H++/3EftTJ\nQkTik8qEkT8PozAZHHQQvP56aC9d6idio3r9e+7xE+APPNC758RZn12zBqZP90n1nnv8KGn1arjr\nLj/H09DQv/1RrTpQLALFInqpTBj5q9UWJowDDoC33grtJUuiSxjg/xK//3649FJ/xFYSLV/ur+Z7\n1FF+dPHcc/7IsWnT/DWQRCRdUpswYPeS1IgR/ozhvCVLdp8Aj8Lkyf7kvksu2fPSxsX6sz67erX/\n3fOTToIjjoCXX4bbb0/OJU5Uqw4Ui0CxiF4qE0b+JLCeRhgrV/rzM6J23XV+gv3uu6N/rb3Ztg3+\n/d/9/MSwYf7w3+uv92fEi4ikMmF0N8LIJwznfN2+u6vXllN9vZ8LuPba3RNWsajrs08/DS0t8Mwz\n/kS7W27xF95LItWqA8UiUCyil+qEUTzCyJekNmzwX+SDB/dPf4491v9i39VX98/rFdqxA777XX8W\n+o03+vMpIvz9FRGpYKk9DwN2Txj5EUZ/ji4K3Xijn2SfO9dfbr1YFPXZt97yh8Lu3OlHFf39nvtK\ntepAsQgUi+ileoRRWJJqbPSJZPPmeBLGkCHwb//mD1/tD8uW+WtaTZoEjz9eOclCROKT6oRRfFLe\n0KH+5zXjSBjgj5Z69VX/BV6snPXZefP8dZ6uuQZuvjmakxOjpFp1oFgEikX0UpkwuipJgf9xnU2b\n/FFLI0b0f7/q6/18wjXX7P0w276aN89fxuPHP4YvfCGa1xCR6pTKhNFVSQr8RfI2bYKNG8OVbPvb\neef561w9VvQr5+Woz778Mnz84z5ZfOYzJe8uNqpVB4pFoFhEL9UJo6sRxubN/iipuM49qKnx15i6\n9dby7nfjRvjEJ/x5FZWcLEQkPqlMGPmSVBJHGOAvvbFihT9yKa+U+qxzcOGF/pLjl11Wev/iplp1\noFgEikX0Upkw9jbCiDth1Nf7S4b85Cfl2d+dd/oLK952W3n2JyLppIRRID/pHWdJKu/ii/2VbDdu\n9O2+1meXLYNvfQvuvbf3v72RdKpVB4pFoFhEL5UJI+klKfC/xnfmmf6Lvq+cgy9+0Z/fEfWFFEWk\n+qUyYfRm0jvuhAHhF+ygb/XZBx7wv11eDfMWhVSrDhSLQLGIXqoTRk8jjLhLUgAf+xi0tfkr5+6r\nzk5/tNVtt1XeiXkikkypTBg9nbj33nv+NmRI//erWH29PwR25sx9r8/ed5+/iOBpp0XStVipVh0o\nFoFiEb1UJozuSlIDB8K6dTBgQHL+Kp82DX796317zo4dcNNNfrJbRKRcUp0wiktSgwaFS5snxSmn\n+DO0H3ww2+vnPPKIn4Op1j+4VKsOFItAsYheKhNGdyWpgQP9/EVxIolTfT2cfjrMmdP75/z0p/Dl\nL4NZdP0SkfRJZcLoriQ1YAC0tydrhAH++k8rV2Z6te2qVfDUU76UVa1Uqw4Ui0CxiF6qE0bxSKKu\nDrZsSdYIA2DqVH8xwh079r7tzJn+AoYDB0bfLxFJl8gShpmdaGZzzGxe7v6EgnXTzazNzJaY2ZSo\n+tCd7kpSdXWwdWvyEsbBB0NjY5bFi/e+7W9/W/0XF1StOlAsAsUielGOMP4DuN45dyxwA/A9ADOb\nAJwPjAfOAu4w699qe3cjjPr6ZCYMgA9+EJ55pudt1q6FJUvgH/6hf/okIukSZcJ4DcifL90ErMk9\n/iQw0zm33Tm3EmgDJkXYjz30NMJIYkkK4FOfyvD00z1v89BD8I//CA0N/dOnuKhWHSgWgWIRvSgT\nxjXA/zGzV/GjjfyvVY8CVhVstya3rN90N+md1JIU+N/f3tsI49FHfcIQEYlCSV+NZjYLGFm4CHDA\n9cDlwOXOud+Z2XnA3cAZ+/oara2tNDc3A9DU1ERLS8uuvyTyNct9bY8b59uLFmWprw/r58/Psm0b\n1NWVtv8o2m++mWX1aj+KOOecPdc7B7NnZ7nwQoD4+xtlO78sKf2Jsz1//nyuuOKKxPQnzvYPfvCD\nsnw/VGI7m80yY8YMgF3fl5FwzkVyA94tam/I3V8DXF2w/BFgcjf7cFFYu9Y5cO4vf9l9+aJFfvnx\nx0fysiWZPXu2O+005/74x67XL1niXHNz//YpLrNnz467C4mhWASKRZD77iz793qUJak2MzsVwMxO\nx89VADwMTDOzBjM7DBgL7MNpaaXr6bDarpYnQSaT4ZhjYOHCrtc/9RR8+MP926e45P/CEsWikGIR\nvSi/Gr8E/NjMGoCtwBcBnHOLzexXwGKgE7g0lxH7TXdzGPkT9pKYMACOOgr+/Oeu16UpYYhIPCIb\nYTjnnnfOTXbOHeuc+5Bzbl7Bupudc2Odc+Odc49F1Yfu1PRwlFThfZJks1mOOqr7EcbTT8PJJ/dv\nn+JSOJeRdopFoFhET2d6F0hywgD/q3nLlsHOnbsv37zZXxJEv6onIlFKdcKopBFGJpNh0CAYNgxW\nr9593Ysv+mSRtGtgRUW16kCxCBSL6KUyYXRXkkr6HAbAuHGwfPnuyxYuhKOPjqc/IpIeqUwYlViS\nytdnx43zP9taaOlSGD++//sUF9WqA8UiUCyil+qEUUklqbzDD4cVK3Zftnw5jB0bT39EJD1SmTDy\nJamuLj7Y1fIkyNdnP/ABWLNm93UrVvhEkhaqVQeKRaBYRC+VCcPM34pHGN2VqpJk9Gh/RFSec+lL\nGCISj1QmDPDJoThh5JNIEhNGvj47evTuR0mtWweNjTBkSDz9ioNq1YFiESgW0Uttwqip6Tox1Ncn\n+/DUUaN8SSp/LsaaNT6JiIhELbUJo6sRBvgkksQRRr4+u99+MHgwvPmmX752rf9FvjRRrTpQLALF\nInqpTRiDBnX9Q0NJTRiF3v9+X4oCnzAOOije/ohIOqQ2YSxe7JNGsfr6ZCaMwvrsiBHw1lv+8Wuv\npW+EoVp1oFgEikX0UpswRozoenkljDCGD4e33/aP01iSEpF4pDZhdCepCaOwPls4wkhjSUq16kCx\nCBSL6ClhFElqSapQ4QjjnXe6Hy2JiJSTEkaRpI4wupvD2LABmpri6VNcVKsOFItAsYieEkaRpCaM\nQiNGhBHG+vXpSxgiEg8ljCJJLUkV1meHD0/3CEO16kCxCBSL6ClhFKmEEcawYT5RdHTA9u3+0iAi\nIlFTwiiS1IRRWJ9tbIT2dti40Y8uzOLrVxxUqw4Ui0CxiJ4SRpGklqQK5RNGGstRIhIfJYwiSR1h\nFNZn054wVKsOFItAsYieEkaRpCaMQmlPGCISDyWMIkktSXU1h7FhAwwdGl+f4qJadaBYBIpF9JQw\nilTCCGPAANi2DTZtgoED4+6NiKSFEkaRpCaMwvqsmR9lbNzY9SXaq51q1YFiESgW0UvgV2O8rroK\nxo2Luxd719joS1JpTBgiEg+NMIqccgoceGDcvdhTcX02zQlDtepAsQgUi+gpYVSoxkZ/Hak0JgwR\niYc55+LuQ7fMzCW5f3E64QT/w0ktLXDjjXH3RkSSxMxwzpX9GhAaYVSoNJekRCQeJSUMMzvPzBaZ\n2Q4zO65o3XQzazOzJWY2pWD5cWa2wMyWmdkPSnn9NNEcRqBadaBYBIpF9EodYSwEPgX8pXChmY0H\nzgfGA2cBd5jtukTeT4AvOOeOAI4wszNL7EMqzJ8/f7d2mhNGcSzSTLEIFIvolZQwnHMvOefagOJa\n2TnATOfcdufcSqANmGRm7wcGO+fm5rb7BXBuKX1Iiw0bNuzWTnPCKI5FmikWgWIRvajmMEYBqwra\na3LLRgGrC5avzi2TfdTYCO+9l86EISLx2OuJe2Y2CxhZuAhwwHXOud9H1THZ3cqVK3dr5380KY0J\nozgWaaZYBIpF9MpyWK2ZzQaucs69kGtfAzjn3K259iPADcArwGzn3Pjc8mnAqc65r3SzXx1TKyLS\nB1EcVlvOS4MUdu5h4D4zuw1fchoLzHHOOTPbaGaTgLnA54EfdrfDKN6wiIj0TamH1Z5rZquAk4A/\nmNmfAJxzi4FfAYuBPwKXFpyBdxnwc2AZ0Oace6SUPoiISP9I9JneIiKSHIk809vMpprZ0tzJfVfH\n3Z+omNlKM/t/ZjbPzObklg0zs8fM7CUze9TMhhZsXzUnQ5rZz83sDTNbULCsbO/dzBrMbGbuOc+Y\n2SH99+72TTexuMHMVpvZC7nb1IJ1VRkLMxttZk+Y2YtmttDMvpZbnrrPRRexuDy3PN7PhXMuUTd8\nElsOHArUA/OBI+PuV0Tv9WVgWNGyW4Fv5h5fDdySezwBmIefd2rOxSg/QnwOODH3+I/AmXG/t168\n9w8DLcCCKN478BXgjtzjC/DnBcX+vvchFjcA3+hi2/HVGgvg/UBL7vEg4CXgyDR+LnqIRayfiySO\nMCbh5zZecc51AjPxJwJWI2PPUd45wH/mHv8n4cTGT1JFJ0M6554C1hctLud7L9zXb4DTy/4myqSb\nWMCeJ8RC306KrYhYOOded87Nzz3eBCwBRpPCz0U3scifsxbb5yKJCaP4pL9qPrnPAbPMbK6Z/c/c\nspHOuTfAf2iA/K9zpOFkyAPL+N53Pcc5twPYYGb7R9f1SHzVzOab2c8KyjCpiIWZNeNHXc9S3v8T\nlRyL53IJlJEWAAAB3klEQVSLYvtcJDFhpMkpzrnjgI8Dl5nZR/BJpFCaj0oo53uvtEO07wDGOOda\ngNeB/13GfSc6FmY2CP8X79dzf11H+X+i0mIR6+ciiQljDVA4+TI6t6zqOOdey92/CfwOX457w8xG\nAuSGk+tym68BPlDw9Hxculteicr53netM7NaYIhz7p3oul5ezrk3Xa64DPxf/GcDqjwWZlaH/4K8\nxzn3UG5xKj8XXcUi7s9FEhPGXGCsmR1qZg3ANPyJgFXFzBpzfz1gZgOBKfir/z4MtOY2uwjI/6d5\nGJiWO7LhMMLJkK8DG81skpkZ/mTIh6gMxp4nfLbmHpf63h/O7QPgn4AnInsX5bFbLHJfjHmfBhbl\nHld7LO4GFjvnbi9YltbPxR6xiP1zEffRAN0cITAVf1RAG3BN3P2J6D0ehj8CbB4+UVyTW74/8Hju\n/T8GNBU8Zzr+6IclwJSC5cfn9tEG3B73e+vl+78fWAt0AK8CFwPDyvXegQH4k0fb8HXw5rjf8z7G\n4hfAgtxn5Hf4On5VxwI4BdhR8P/ihdx3Qdn+T1RBLGL9XOjEPRER6ZUklqRERCSBlDBERKRXlDBE\nRKRXlDBERKRXlDBERKRXlDBERKRXlDBERKRXlDBERKRX/j/BdqLLYWkp5AAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x106a911d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "bb, aa = signal.ellip(6, 1, 60, 1300.0/24000, 'low')\n",
    "w,h = signal.freqz(bb, aa)\n",
    "w = w * 48000 / (2 *np.pi)\n",
    "pl.ylim((-100, 10))\n",
    "pl.plot(w, 20 * np.log10(abs(h)));\n",
    "bb,aa"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYMAAAD7CAYAAACIYvgKAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XmQVPW99/H3F9lBGFwAHZYBWRyQZDSAxiU2EAjqFair\nlZikbjKJGnPJ9YlGc5XHVHhukvuo2Z5YKXPLysa9icaosQwqGnHplEsAF0ZEwBll35RFURZR4Pf8\n8etxGhxmpqf7zDnn159X1anpc7p7zvc7P+hvn9+3+xxzziEiIuWtU9wBiIhI/FQMRERExUBERFQM\nREQEFQMREUHFQEREgM5xB2Bm+myriEg7OOesVL8rEUcGzrlgl7lz58Yeg/JTfuWYX8i5OVf699CJ\nKAYhW7t2bdwhREr5pVvI+YWcWxRUDERERMUgarW1tXGHECnll24h5xdyblGwKOaeCgrAzMUdg4hI\n2pgZLrQGcsiy2WzcIURK+aVbyPmFnFsUVAxERETTRCIiaaRpIhERKTkVg4iFPm+p/NIt5PxCzi0K\nkRcDM5tuZqvMrN7Mboh6fyIiUrhIewZm1gmoB6YAm4Hngcucc6vyHqOegYhIgdLWM5gINDjn1jnn\nPgTuBmZGvE8RESlQ1MWgEtiQt74xt61shD5vqfzSLeT8Qs4tCrGfwjoOzsGuXbBxI7z5Jrz11uHL\nzp3w3nuHL/v2+ecdOuR/Ogdm0KuXX3r3blpOOAEGDPDLjh1w4AAMGgRDhkDPnnFnX5gPPoDVq+G1\n15qW+nrYssXfv2+fzx+ge3cYNQrGjoUxY/zPUaOga9fSxLJ9O9xxB9x7L9x1l9+HiJRG1MVgEzAk\nb31QbtthamtrqaqqAqCiooKamhoymQzQVN3buv7EE1l27YKhQzNs2gRPPpll2zY45hi/Xl/v17t2\nzVBZCd27Z+nXD047LcOAAdCnT5ZBg+CsszIceyysWpWlZ0+YNClDp07w3HN+f+edl+HQIb+/99+H\nMWMy7N7t73/nHejbN8O6dfDKK/Dww1n27MmwYYPf34ABfn9Dh8IHH2QZOBAuusivv/xyFrO251uK\ndeegujrDa6/B/PlZNmyAvXsz1NfD2rVZTjwRTj89w6hR0Ldvlksu8fECLFoEkGXixAz79sH992dZ\nswaWL8/w6quwenWWUaNg1qwMkybBgQNZunUrLL41a+C55zLcdx98+tNZPvUpuPTSDEuWwAsvRPv3\nadzWkePRkesh55fJZBIVT7Hr2WyWefPmAXz0ellKUTeQjwFewzeQtwBLgC8651bmPabFBvLBg03v\nzrdv9+/k85e33jp8fccOqKjw78QrK/3PxqVxvbIS+vSJLO2jOnTIx7h+Paxb17Tkrx886I8ghg5t\nWoYMgf794cQT/VHHCSdAt26F7XfHDv9ufvNmv2zc6N/hN77T79oVRo/27+RHj25aTjmluHf2+/bB\nc8/BU0/Bk0/CsmUwfjxMmODf2Y8ZA9XVTePxwQfQ0ACrVsHKlfD3v8Py5TB7Nlx1lf87AHz96/6x\nf/iDP0ITKTelbiBH/g1kM5sO3IbvT/zWOXfLEfe78eMdBw7Ahx/y0c/du30B2L/fT0Mce6x/Eezf\nv2kKprnlxBOhS5dIUypI/ruutti1q/lisW2bX7Zv90v37v7v0q2bv92tW1OB2Lfv8GX3bv/3O/lk\nOOkk/7OyEkaObHrRP+64jsnvvffgmWegrg5WrPAv+KtWQb9+0KOHz3XoUDj1VL+ccQbMmvXx4rd3\nL5x1FnzrW75IRKXQ/NIm5PxCzg1KXwwi7xk45x4FRrf0mNtvh86d/Yt448/evf0LWM+e5fXOr29f\nGDfOL0fjHLz7LuzZ44vl++/7n/v3+/t79PBL9+7+Z+/e/nYSHHssXHCBXxodOuSLwN69/kikLUc9\nPXvCfffBOef4o4wzzoguZpFyoHMTSardcw/MmQMvvuinB0XKReqmiVoNQMVAinT11b4Hcv/95XUU\nKeUtbV86K3uNnwYIVRLy++lPYc0aePDB0v/uJOQXpZDzCzm3KKgYSOp16wZf+hIsXBh3JCLppWki\nCcKSJXD55f57HSLlQD0DkWYcOADHHw+vv+4/XiwSOvUMUib0ecuk5Ne5M5x7rv+SWiklJb+ohJxf\nyLlFQcVAgjFpkv+ms4gUTtNEEowXXoCvfhVefTXuSESip56ByFEcPOj7BvX1TecwEgmVegYpE/q8\nZZLyO+YYOO88KGVIScovCiHnF3JuUVAxkKBMmlTaYiBSLjRNJEF56SX48pf92VBFQqZpIpEWfPKT\nsHWrX0Sk7VQMIhb6vGXS8jvmGPjMZ0o3VZS0/Eot5PxCzi0KKgYSHPUNRAqnnoEEp64OvvAFf0lP\nkVCpZyDSik98wl8adPPmuCMRSQ8Vg4iFPm+ZxPw6dSpd3yCJ+ZVSyPmFnFsUVAwkSOobiBRGPQMJ\n0rJlcMkl0NAQdyQi0VDPQKQNTjsN3n7bXxtZRFqnYhCx0Octk5pfp05w/vnFTxUlNb9SCTm/kHOL\nQlHFwMwuNbPlZnbQzM444r45ZtZgZivNbFpxYYoUTn0DkbYrqmdgZqOBQ8AdwPXOuZdy26uBu4AJ\nwCDgcWBkc80B9QwkKsuXw8yZ8MYbcUciUnqJ6hk4515zzjUARwY0E7jbOXfAObcWaAAmFrMvkUKN\nHQvvvQfr18cdiUjyRdUzqAQ25K1vym0rO6HPWyY5P7Pi+wZJzq8UQs4v5Nyi0Lm1B5jZQmBA/ibA\nATc55x4sRRC1tbVUVVUBUFFRQU1NDZlMBmga0LSu19XVJSqecsuvsjLLn/4EX/lKmPmFPn5ab1rP\nZrPMmzcP4KPXy1IqyfcMzOwp4Lq8nsGNgHPO3ZpbfxSY65xb3Mxz1TOQyKxYARddBGvWxB2JSGkl\nqmdwhPyg5gOXmVlXMxsGjACWlHBfIm1SXQ1798LatXFHIpJsxX60dJaZbQDOAh4ys0cAnHMrgHuA\nFcACYHa5vv1vPMwLVdLzM4NMpv19g6TnV6yQ8ws5tygU+2miB5xzg51zPZxzJznnLsi772bn3Ajn\nXLVz7rHiQxVpH33fQKR1OjeRBG/VKvjc5/xUkZVshlUkXknuGYgk0ujR8MEH6huItETFIGKhz1um\nIb/GvsFTTxX+3DTkV4yQ8ws5tyioGEhZUN9ApGXqGUhZqK+HKVP8qSnUN5AQqGcg0g4jR8KhQ7B6\nddyRiCSTikHEQp+3TEt+7e0bpCW/9go5v5Bzi4KKgZQN9Q1Ejk49Aykbr7/uz2K6caP6BpJ+6hmI\ntNMpp/jLYb7+etyRiCSPikHEQp+3TFN+7ekbpCm/9gg5v5Bzi4KKgZSVSZPa9+UzkdCpZyBlZd06\nmDABtm71U0YiaaWegUgRhg6F44+HpUvjjkQkWVQMIhb6vGUa85s+HR59tG2PTWN+hQg5v5Bzi4KK\ngZSdQoqBSLlQz0DKzr590L8/bNgAFRVxRyPSPuoZiBSpRw8491x44om4IxFJDhWDiIU+b5nW/No6\nVZTW/Noq5PxCzi0KKgZSlqZPh7/9DTRDKeKpZyBlyTkYPhwefhjGjIk7GpHCqWcgUgJm+lSRSD4V\ng4iFPm+Z5vzaUgzSnF9bhJxfyLlFoahiYGY/NrOVZlZnZn8xsz55980xs4bc/dOKD1WktCZPhn/8\nA/bsiTsSkfgV1TMws88CTzrnDpnZLYBzzs0xszHAncAEYBDwODCyueaAegYSp0mT4LvfhQsvjDsS\nkcIkqmfgnHvcOXcot7oI/8IPMAO42zl3wDm3FmgAJhazL5EoqG8g4pWyZ/B1YEHudiWwIe++Tblt\nZSf0ecu059daMUh7fq0JOb+Qc4tC59YeYGYLgQH5mwAH3OScezD3mJuAD51zf2pPELW1tVRVVQFQ\nUVFBTU0NmUwGaBrQtK7X1dUlKh7ld/j6zp1Ztm+HN97IcMop4eUX+viV03o2m2XevHkAH71ellLR\n3zMws1rgSmCyc25/btuN+P7Brbn1R4G5zrnFzTxfPQOJVW0tTJwIs2fHHYlI2yWqZ2Bm04HvAjMa\nC0HOfOAyM+tqZsOAEcCSYvYlEhX1DUSK7xn8EugNLDSzl8zsVwDOuRXAPcAKfB9hdrm+/W88zAtV\nCPlNnQp//zu8//7H7wshv5aEnF/IuUWh1Z5BS5xzI1u472bg5mJ+v0hHOP54GDcOsll/lCBSjnRu\nIhHgxz/210e+/fa4IxFpm1L3DFQMRIBVq/x00fr1/rxFIkmXqAaytC70ectQ8hs9Grp3h9wnLT8S\nSn5HE3J+IecWBRUDEfzRwIwZ8OCDcUciEg9NE4nkZLNw/fXwwgtxRyLSOvUMRCLy4YcwcCAsWwaV\nZXnyFEkT9QxSJvR5y5Dy69IFLrgAHnqoaVtI+TUn5PxCzi0KKgYieS6+GObPjzsKkY6naSKRPLt2\nweDBsGUL9OoVdzQiR6dpIpEI9e0LZ54JCxfGHYlIx1IxiFjo85Yh5pc/VRRifvlCzi/k3KKgYiBy\nhIsv9k3kgwfjjkSk46hnINKMcePgjjvg7LPjjkSkeeoZiHQAfRtZyo2KQcRCn7cMNb8ZM3zfINT8\nGoWcX8i5RUHFQKQZEybAjh2waVPckYh0DPUMRI7iyiuhuhq+8524IxH5OPUMRDrIrFnwwANxRyHS\nMVQMIhb6vGXI+U2ZAi++mOWtt+KOJDohj1/IuUVBxUDkKLp3h/HjDz9xnUio1DMQacGdd8Kf/6yT\n10ny6HoGIh3o7bdh6FDYvBl69447GpEmaiCnTOjzlqHn9/LLWc48Ex57LO5IohHy+IWcWxSKKgZm\n9gMze9nM6szscTMblHffHDNrMLOVZjat+FBF4qFPFUk5KGqayMx6O+d2525fDXzCOXelmY0B7gQm\nAIOAx4GRzc0HaZpIkm7DBqipga1b/dXQRJIgUdNEjYUgpxewI3d7BnC3c+6Ac24t0ABMLGZfInEZ\nPBiGD4dnnok7EpHoFN0zMLMfmdl6oBa4Obe5EtiQ97BNuW1lJ/R5y3LJb+bMMKeKQh6/kHOLQufW\nHmBmC4EB+ZsAB9zknHvQOfc94HtmdgPwC+BrhQZRW1tLVVUVABUVFdTU1JDJZICmAU3rel1dXaLi\nUX7ty2/WrAwXXQSzZmUxS058Gr/yWc9ms8ybNw/go9fLUirZR0vNbDCwwDk3zsxuBJxz7tbcfY8C\nc51zi5t5nnoGknjOwciRcO+9cPrpcUcjkrCegZmNyFudBdTlbs8HLjOzrmY2DBgBLClmXyJxMtOn\niiRsxfYMbjGzZWa2FMgA1wE451YA9wArgAXA7HJ9+994mBeqcsovxGIQ8viFnFsUWu0ZtMQ5d2kL\n991MU0NZJPU+/WnYsgVWr/afLhIJiU5HIVKAK66AsWPh2mvjjkTKXaJ6BiLlJsSpIhFQMYhc6POW\n5ZbflClQVwfbtsUTT6mFPH4h5xYFFQORAvToAVOn6hoHEh71DEQK9Mc/+u8b/PWvcUci5UzXMxCJ\nWeM1DrZsgV694o5GypUayCkT+rxlOebXrx9MnBjGNQ5CHr+Qc4uCioFIO+hTRRIaTROJtEPjNQ7e\nfBM6F/XVTZH20TSRSAIMHgxVVfD003FHIlIaKgYRC33espzzmzUr/Z8oCnn8Qs4tCioGIu3U2DfQ\nLKeEQD0DkXZyDkaMgL/8xfcPRDqSegYiCaFrHEhIVAwiFvq8Zbnnl/ZiEPL4hZxbFFQMRIpw9tmw\neTOsWRN3JCLFUc9ApEiXXw7jxsE118QdiZQT9QxEEmbmzPR/xFRExSBioc9bKj9/jYMXXoBdu6KP\np9RCHr+Qc4uCioFIkXr1gnPOCePEdVK+1DMQKYFf/hJeegl+//u4I5FyoesZiCTQ6tVNnyzqpONt\n6QBqIKdM6POWys8bPtxf5+DFF6ONp9RCHr+Qc4tCSYqBmV1nZofM7Li8bXPMrMHMVprZtFLsRyTJ\nLroIHn447ihE2qfoaSIzGwT8BhgNfMo5t9PMqoG7gAnAIOBxYGRz80GaJpJQPPkk3HgjLFkSdyRS\nDpI4TfT/gO8esW0mcLdz7oBzbi3QAEwswb5EEuvcc6G+3l/wRiRtiioGZjYD2OCce+WIuyqBDXnr\nm3Lbyk7o85bKr0nXrvDZz8Ijj0QXT6mFPH4h5xaFVi/YZ2YLgQH5mwAHfA/438DUYoOora2lqqoK\ngIqKCmpqashkMkDTgKZ1va6uLlHxKL9o8zvllCy//z3U1iYj/nIfv5DWs9ks8+bNA/jo9bKU2t0z\nMLPT8L2AvfgCMQh/BDAR+DqAc+6W3GMfBeY65xY383vUM5BgbN0K1dXw1lvQpUvc0UjIEtMzcM4t\nd84NdM4Nd84NAzYCpzvn3gLmA18ws65mNgwYAaitJsEbONBf8ObZZ+OORKQwpfyegcMfIeCcWwHc\nA6wAFgCzy/Xtf+NhXqiU38ddeGF6PmIa8viFnFsUSlYMckcIO/PWb3bOjXDOVTvndNYWKRv6voGk\nkU5HIVJihw7BSSfBokUwbFjc0UioEtMzEJHmdeoEF1wACxbEHYlI26kYRCz0eUvl17y0TBWFPH4h\n5xYFFQORCEydCk8/Dfv3xx2JSNuoZyASkYkT4Sc/gfPPjzsSCZF6BiIpMWkSPPVU3FGItI2KQcRC\nn7dUfkeXhmIQ8viFnFsUVAxEInLuuf5iN/v2xR2JSOvUMxCJ0Nlnww9/CFOmxB2JhEY9A5EUScNU\nkQioGEQu9HlL5deypBeDkMcv5NyioGIgEqGzz4aXX4bdu+OORKRl6hmIROwzn4GbboLPfS7uSCQk\n6hmIpEzSp4pEQMUgcqHPWyq/1iW5GIQ8fiHnFgUVA5GInXUWvPoqvPtu3JGIHJ16BiIdYPJkuO46\nfzZTkVJQz0AkhZI8VSQCKgaRC33eUvm1TVKLQcjjF3JuUVAxEOkAEydCfT28/XbckYg0Tz0DkQ4y\nbRrMng2zZsUdiYRAPQORlErqVJEIqBhELvR5S+XXdkksBiGPX8i5RaGoYmBmc81so5m9lFum5903\nx8wazGylmU0rPlSRdBs/Htatg+3b445E5OOK6hmY2VzgPefcz4/YXg3cBUwABgGPAyObaw6oZyDl\nZPJkuOEGnadIipfEnkFzwcwE7nbOHXDOrQUagIkl2JdIqtXUQF1d3FGIfFwpisG/mVmdmf3GzPrm\ntlUCG/Iesym3reyEPm+p/AqTtGIQ8viFnFsUOrf2ADNbCAzI3wQ44CbgV8APnHPOzH4E/Ay4otAg\namtrqaqqAqCiooKamhoymQzQNKBpXa/L/c9PSjzKL978DhzI8uyzAGHmp/Xo1rPZLPPmzQP46PWy\nlEr2PQMzGwo86Jz7hJndCDjn3K25+x4F5jrnFjfzPPUMpGx8+CH07QvbtkGvXnFHI2mWqJ6BmQ3M\nW/1nYHnu9nzgMjPrambDgBHAkmL2JRKCLl2guhpeeSXuSEQOV2zP4MdmtszM6oDzgWsBnHMrgHuA\nFcACYHa5vv1vPMwLlfIrXJL6BiGPX8i5RaHVnkFLnHNfaeG+m4Gbi/n9IiFKUjEQaaRzE4l0sGee\ngeuvh0WL4o5E0qzUPQMVA5EO9u67cPLJsGsXHHNM3NFIWiWqgSytC33eUvkVrk8fGDjQn9I6biGP\nX8i5RUHFQCQG6htI0miaSCQGP/oRvPce3Hpr3JFIWmmaSCQAOjKQpFExiFjo85bKr31OPx2WLoW4\nD4pDHr+Qc4uCioFIDE4+2ReCrVvjjkTEU89AJCZTp8K118KFF8YdiaSRegYigVDfQJJExSBioc9b\nKr/2S0IxCHn8Qs4tCioGIjE5/fT4i4FII/UMRGJy4IC/tsHWrXDssXFHI2mjnoFIIDp3hrFjYdmy\nuCMRUTGIXOjzlsqvOHH3DUIev5Bzi4KKgUiM4i4GIo3UMxCJ0XPPwbe/Dc8/H3ckkja6noFIQHbv\nhgED/LUNOhd13UEpN2ogp0zo85bKrzi9e8OgQfDaa5Hu5qhCHr+Qc4uCioFIzD75SX2iSOKnaSKR\nmF19NYwa5X+KtJWmiUQC06ePvy6ySJyKLgZmdrWZrTSzV8zslrztc8ysIXfftGL3k1ahz1sqv+LF\nWQxCHr+Qc4tCUZ9fMLMMcDEwzjl3wMxOyG2vBj4PVAODgMfNbKTmg0Q+rk8fWLs27iik3BXVMzCz\nPwN3OOeePGL7jYBzzt2aW38E+D/OucXN/A7VCClrd94JCxb4nyJtlbSewSjgM2a2yMyeMrNP5bZX\nAhvyHrcpt01EjqCegSRBq9NEZrYQGJC/CXDA93LP7+ecO8vMJgD3AsOjCDStstksmUwm7jAio/yK\nV+pisG4d3HEH9OwJU6bA+PHQpUvzjw15/JKcm3PwzjvNL3v2wL59sHev/7lvnz/DrXNNSxRaLQbO\nualHu8/Mvgncn3vc82Z20MyOxx8JDMl76KDctmbV1tZSVVUFQEVFBTU1NR8NYmMTKK3rdbkTzyQl\nHuWXvPwaGuDdd4v/fc8/DzfckOXFF+HKKzPs3w//8i9Ztmzx90+eDIMHZ+nfv3zGL45152Ds2Axr\n1sBDD/m/f+fOGbZsgfr6LDt2wK5dGbp1g549s/TuDYMHZ6iogH37svToASNGZOjRA3bsyNKtG4wZ\nk6G+PsvixfMAOOGEKkqt2J7BN4BK59xcMxsFLHTODTWzMcCdwJn46aGFQLMNZPUMpNy98QZMm+Z/\ntsff/gb/+Z+wfj1ccw1cfvnh10fYvh2yWVi4EO6/H267Db70pZKEXtacgy1b4JVXDl/q66FbNxg2\nDIYP98uQIXDSSX4ZONCfgqRHj+L2n6hzE5lZF+B3QA2wH7jOOff33H1zgMuBD4FvO+ceO8rvUDGQ\nsrZtG1RX+xftQq1eDRMnwu23wyWXtH5+o6VL4YtfbHqOLqrTdu+/Dy++6E8u+Oyz/qdzMG7c4cvo\n0VBREX08iSoGJQkg8GKQTfC8ZSkov+Lt3+9flPfvByvwv/ZVV/l3mT/4Qdufs2ePP4J46im47ros\n//qvmcJ2mhLFjp1z/jQhDz0EjzziC+mpp8I55/jl7LP9eaUKHbNSKXUx0HkSRWLWrZt/Qdm/H7p3\nb/vzNm2C++4r/CR3vXrBr3/tn3vFFX6/11xT2O8I1f79fjrtoYfg4Yeha1e4+GL4/vf9i3/v3nFH\nGB0dGYgkwIknwquvQv/+bX/Od77ji8jPftb+/a5fD5Mnw7//O3zjG+3/PWlXVwe/+x386U9+ym7G\nDPinf/JTPnG982+NjgxEAtT48dK2FoPt22HePFi+vLj9DhkCjz4K550HJ5/sXwDLxdtvwx//6IvA\nzp1QWwtLlvjGbznSieoi1vjRs1Apv9Io9LsGt90Gn/+8fwEvRjabZcQIeOAB+NrX/IthKI42dtu3\nw003wciR8I9/wE9+AmvWwH/8R/kWAtCRgUgiFFIMdu2C//qv0r5wn3mmf4c8cyY8/TSMGFG6350U\nb77pp9R+8xtfSF94AXJfbxLUMxBJhIsvhiuv9HPVrbnlFt9f+MMfSh/HHXf4d8rPPVdY/yLJ9uzx\nn7b69a/hy1/2/ZHBg+OOqnjqGYgEqK1HBnv3wi9+AU88EU0cV10FGzb44vT00/7TNGn22GPwzW/6\nTwItX178tFrI1DOImObU0y1pPYPf/ta/sI0dW5r9NpffD38Ixx8PP/95afYRh+3bYdq0LN/4Bvzq\nV75RrELQMhWDiDWe+yVUyq802lIMPvjAT+HMmVO6/TaXn5n/dvJPf+q/4Zw2d98Np50G+/fXsXw5\nTJ8ed0TpoGIQsXfeeSfuECKl/EqjLcVgxQr/pacJE0q336PlN2wYXH89fOtb0Z0ls9QOHvQxf//7\nMH8+TJr0TtBfEis1FQORBGhLMdi50596oqNcd53vH9x7b8fts7327IFLL/XnDlq0yJ97SQqjYhCx\ntYFfz1D5lUafPv4joy3ZuRP69SvtflvKr0sX/+mia69tPbY4bd4M558Pffv6M7ged5zfHvq/zVJL\nxEdLYw1ARCSlgjprqYiIxE/TRCIiomIgIiIxFwMzm25mq8ys3sxuiDOW9jKztWb2spktNbMluW39\nzOwxM3vNzP5mZn3zHj/HzBrMbKWZTYsv8uaZ2W/N7E0zW5a3reB8zOwMM1uWG9tfdHQeR3OU/Oaa\n2UYzeym3TM+7L235DTKzJ83sVTN7xcz+V2576sewmdyuzm0PYvzMrJuZLc69lrxqZv83t71jxs45\nF8uCL0SvA0OBLkAdcGpc8RSRx2qg3xHbbgX+PXf7BuCW3O0xwFL8aUCqcvlb3DkcEfu5+MuYLism\nH2AxMCF3ewHwubhzayG/ucB3mnlsdQrzGwjU5G73Bl4DTg1hDFvILaTx65n7eQywCDino8YuziOD\niUCDc26dc+5D4G5gZozxtJfx8SOsmcB/527/NzArd3sGcLdz7oBzbi3QgP87JIZz7hng7SM2F5SP\nmQ0EjnXOPZ973P/kPSdWR8kP/DgeaSbpy2+rc64ud3s3sBIYRABjeJTcKnN3hzJ+e3M3u+FfV96m\ng8YuzmJQCWzIW99I08CmiQMWmtnzZnZFbtsA59yb4P8BA43nfzwy502kI+f+BeZTiR/PRmkY238z\nszoz+03eYXiq8zOzKvxR0CIK/zeZ6Bzzcluc2xTE+JlZJzNbCmwFss65FXTQ2KmBXLxznHNnABcC\n3zKz8/AFIl9on98NLZ9fAcOdczX4/4RFXEgyGcysN3Af8O3cu+hg/k02k1sw4+ecO+ScOx1/NHee\nmWXooLGLsxhsAobkrQ/KbUsV59yW3M9twAP4aZ83zWwAQO6Q7a3cwzcB+WdST0vOheaTqjydc9tc\nbnIV+DVNU3epzM/MOuNfLP/gnPtrbnMQY9hcbqGNH4Bz7l38XP94Omjs4iwGzwMjzGyomXUFLgPm\nxxhPwcysZ+5dCmbWC5gGvILPozb3sK8Cjf8h5wOXmVlXMxsGjACSeKFB4/A52ILyyR3K7jKziWZm\nwFfynpMEh+WX+w/W6J+BxisLpzW/3wErnHO35W0LZQw/llso42dmJzROcZlZD2AqvkHcMWMXc+d8\nOv4TAQ0xokWvAAAApUlEQVTAjXHG0s74h+E/BbUUXwRuzG0/Dng8l9tjQEXec+bgu/4rgWlx59BM\nTncBm4H9wHrga0C/QvMBPpX7mzQAt8WdVyv5/Q+wLDeWD+DnaNOa3znAwbx/ly/l/p8V/G8yaTm2\nkFsQ4weMy+W0FHgZuD63vUPGTqejEBERNZBFRETFQEREUDEQERFUDEREBBUDERFBxUBERFAxEBER\nVAxERAT4/zobbTYV19epAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x106b91d50>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "pl.xlim(0, 3000)\n",
    "pl.ylim((-65, 5))\n",
    "pl.plot(w, 20 * np.log10(abs(h)));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "### 双二次(biquad)に分解\n",
    "\n",
    "特性を確認できたら、つづいて実装のための係数を求める。\n",
    "\n",
    "実装では2次フィルタを直列にカスケードする。そのため6次フィルタを3つの2次フィルタに分解することを考える。具体的にはba形式(多項式の係数)で表現された伝達関数(6次式)を、3つの2次式の積に因数分解することになる。しかし、多項式に展開済みの伝達関数の分解は容易ではないので、分解が容易になるようzpk形式という零点と極の位置とゲインで表現した形式で再設計してみる。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([ 0.75834396+0.65185462j,  0.95713883+0.28962951j,\n",
       "         0.97287413+0.23133508j,  0.75834396-0.65185462j,\n",
       "         0.95713883-0.28962951j,  0.97287413-0.23133508j]),\n",
       " array([ 0.95327800-0.05347031j,  0.96605270-0.13301552j,\n",
       "         0.97881481-0.16771403j,  0.95327800+0.05347031j,\n",
       "         0.96605270+0.13301552j,  0.97881481+0.16771403j]),\n",
       " 0.0010766926489612797)"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "z, p, k = signal.ellip(6, 1, 60, 1300.0/24000, 'low', output='zpk')\n",
    "z,p,k"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "ちなみにzpk2tf関数を使うと、zpk形式を伝達関数の係数形式に変換できる。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(array([ 0.00107669, -0.00578907,  0.01354389, -0.0176606 ,  0.01354389,\n",
       "        -0.00578907,  0.00107669]),\n",
       " array([  1.        ,  -5.79629102,  14.04709981, -18.2174942 ,\n",
       "         13.33384205,  -5.22208148,   0.85492756]))"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "signal.zpk2tf(z, p, k)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "### 順序を決定\n",
    "\n",
    "2次フィルタを構成するために、共役な零点zのペアと極pのペアをそれぞれ選択する。共役なzとpは3個毎に並んでいるようなので、それらをペアにして、それぞれ個別に組み合わせて伝達関数を計算する。それらを重ねてプロットし、検討する。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAEACAYAAACgS0HpAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3XeYVdW9//H3dzoMHQaQbkNBUCCKKBbsxthiu0FvFE30\nRk2uRK8/NCFXE2NBY4lGfDTGYIyKijeKsWIZNWgUhRGkiWVoUqUztGHW7491zuwzlSlnz2mf1/Ps\nZ3Y7+6zzZTjfWWvttbY55xAREdmTrEQXQEREUoMShoiINIgShoiINIgShoiINIgShoiINIgShoiI\nNEioCcPMepnZ22Y218zmmNl/R/Z3NLM3zGyhmb1uZu3DLIeIiDSfhTkOw8y6A92dcyVm1gb4FDgL\nuBT4zjl3p5mNAzo6524IrSAiItJsodYwnHMrnXMlkfUtwHygFz5pPB457XHg7DDLISIizRdqDaPK\nG5n1A4qBQcBS51zHmGPrnHOdWqQgIiLSJC3S6R1pjpoCXBOpaVTPUpqfREQkyeWE/QZmloNPFk84\n516M7F5lZt2cc6si/Ryr63itEomISBM45yze12yJGsZjwDzn3B9j9k0FxkTWLwFerP6iKOecFue4\n6aabEl6GZFkUC8VCsah/CUuoNQwzGwlcBMwxs1n4pqdfAROAZ83sMmAxcEGY5UgHpaWliS5C0lAs\nAopFQLEIX6gJwzk3Hciu4/CJYb63iIjEl0Z6p4gxY8YkughJQ7EIKBYBxSJ8LXZbbVOYmUvm8omI\nJCMzw6Vop7fEQXFxcaKLkDQUi4BiEVAswqeEISIiDaImKRGRNKMmqabYvRt27Up0KURE0kJ6J4zO\nneGQQxJdirhQ+2xAsQgoFgHFInyhTw2SMBUVsHEjtGqV6JKIiKSF9O3DWLUKevXyCWPTpvgWTEQk\niakPo7GWLoVBg2DbNtixI9GlERFJeemdMPr2hS5dYO3aRJem2dQ+G1AsAopFQLEIX/omjCVLoHfv\ntEkYIiKJlr59GGPH+oTx8svwq1/BiZrrUEQyg/owGuurr2C//VTDEBGJk/RNGF9+CfvuC0VFsGZN\nokvTbGqfDSgWAcUioFiELz0Txu7dUFoK++yjGoaISJykZx9GaSmMHAnLl8Of/gTz5sHEiXEvn4hI\nMlIfRmOUlARTgqiGISISF+mbMIYO9evqw0g7ikVAsQgoFuFLm4ThnOM3b/+GM54+g12ffgxDhvgD\nRUWqYYiIxEHa9GE88/kz3PLeLQzrNoSHLn2ewvlfQs+e8O23MGwYrFwZcmlFRJJDWH0YaZEwyivK\nGfjgQCb+YCIDl++k7Jwz6btyG7nZubBzJxQW+p8W9/iJiCQddXrX4/GSx+nVrhcn7H0CPabPpmRQ\nZ95d/K4/mJcHrVvDhg2JLWQzqX02oFgEFIuAYhG+lE8YO8p38Nt3f8utx9+KmcGUKWw6/STeLX03\nOEn9GCIizZbyCePhTx/mkO6HcETvI2DOHFixgt6nX8R7S94LTurSJeXvlBo1alSii5A0FIuAYhFQ\nLMKX0k/cW1u2llvfv5VpP57mdzzwAFxxBSP6HcWnz3/K9vLtFOQUqIYhIhIHKV3DGDdtHKMHjebg\nbgf7yQaffx6uvpq2+W0ZWDSQGctn+BPToIah9tmAYhFQLAKKRfhSNmG8sOAF3vrmLX533O/887t/\n+lO48UafHIBj+h4TdHyrhiEi0mwJSxhmdqqZLTCzL8xsXGNeO3f1XK546QqeO/852uW3g7vu8rfN\n/vKXlecc1ecopi+d7jfSoIah9tmAYhFQLAKKRfgSkjDMLAv4E3AKcBAw2swObMhrZ6+azUlPnMR9\np97HYT0Pgz//GR56CCZPhuzsyvOO6HUEHy79kApXkTbTg4iIJFKiahjDgUXOucXOuV3AZOCs+l5Q\n4Sp4dOajnPC3E7j3lHu5cN+z4ZprYMIEeOMN/3S9GN3adKOosIh5a+ZB166wenV4n6YFqH02oFgE\nFIuAYhG+RN0l1RNYGrO9DJ9Eali+aTmvfvkqD854kLzsPIrPeYmDXp8Jpw+EESNgxgzo2LHWNxnZ\neyTTl0xn0F6HwYoV8f8UIiIZJCFTg5jZucApzrkrItv/CQx3zv13tfPcITd25KzcwZzvBnDQ56ux\n4mI47jjfX3HUUfW+zyOfPsK/lvyLv42Y4CcjXLUqrI8kIpI0wpoaJFE1jOVAn5jtXpF9NQy8azvr\nC77hj6030/17Aznh0UcZdc45QFAFjXZ2Vd/OXZLLm2+/CWd2hXXrKH7rLcjOrvN8bWtb29pOxe3i\n4mImTZoEQL9+/QhLomoY2cBC4ARgBfAxMNo5N7/aee7Xv3Z09d/3PPww/OY3cNVVDXufCldB5zs7\nM//q+XTffyh88omfwTYFFRcXV/6iZDrFIqBYBBSLQFrVMJxzu83s58Ab+I73v1RPFlG//32wfskl\ncOSRcPDBe2yNAiDLsjii1xF8sPQDztlrL9+PkaIJQ0Qk0VJuevPJk+Gee+Cjjxo2W/mt793Kum3r\nuPv+BfCzn8EZZ4RUWhGR5KDpzSMuuAB27ICXXmrY+SP7jPQD+KI1DBERaZKUSxhZWXDzzXDrrQ07\nf3jP4cxZPYddXbukdMKIdnCJYhFLsQgoFuFLuYQBcOaZ/smrn32253Nb57ZmYNFAFrfemdIJQ0Qk\n0VKuDyPq5pv9fIJ/+tOerzP2tbEcM/M7zvloE7z4YnwLKSKSZNSHUc1ll8HTT0NZ2Z7PHdl7JNN3\nfa0ahohIM6RswujTBw4/HKZM2fO5R/c9mmnbPselcMJQ+2xAsQgoFgHFInwp/cS9yy+He++Fiy+u\n/7zubbpT0a0rblUpVlHhe84zgXNQXu6XXbuC9epL9WO7d/tnjMQuztXc15hlT693ruoSLX/1defg\n66/hgw8adm5Dj4V5nbr+beKx/9tv4e9/b9xrmvI+TXlNS7//qlXwyCOJe/94y8724wiSSMr2YYB/\nBEavXv67Y7/96r/W5VMv5/4fP0WrL0v9dOeJtHOn74BZtw42b4ZNm/zP6suWLbB9u7+PePv2qktd\n+2K/+CsqICcnWHJzq27Xdiw72//MygoWs6rbjV329ProcfDr0SV2u75jjTk3UdepTV3HGrs/mV+T\n6e/fHGZw3nlNfGkajfSOl7w8uOgimDSp6ojw2hzT9xjWtH+WPitWhJcwnINly6C0FBYvDpYVK/zz\nONau9T/LyvxDnTp1gnbtoG3bmkuHDn5UeqtWUFAA+fn+Z3SJ3Y6u5+X5L/7ol3/0y1hEJA5SuoYB\nMHs2/OAH/js65vlJNSzesJjS4ftzzP1TsVNPbX7hnIP58+G992DWLJgzBz7/HFq3hn32gb59g6VH\nD5+kiop8oujQodFf5JonJ6BYBBSLgGIRUA2jDgcf7J+P9NZbcPLJdZ/Xt0NfPmmfy8pFs9irqQlj\n+3Z45RV47jl4+20oLIRjj4VDD4ULL4RBg6Bz56ZdW0QkyaV8DQP8WIzp0/1ttvV5+YcHUdTzAIb/\n6f8aV5AlS+APf/Cdi0OGwOjRcMop/lYtEZEko3EY9bjwQnj1VVi/vv7zOu19EBu+mdfwC5eVwfXX\nw9Chvo9g9mxfs7j8ciULEck4aZEwOnXyzVF7ugOt74AR7Fy+tP6TombOhMGD/W2L8+fDnXf6W7IS\nRPeYBxSLgGIRUCzClxYJA/zI70cfrf+W6L32H0aXTbtYvGFx/Rd74w3f5HTHHfDkk76TREQkw6VF\nHwb4sWb9+/vv9xEj6jjpiy9Yecwwpr3+ED8+5Me1nzNjBpx2GvzjHw17SpOISJJRH8YeZGfD1VfD\nAw/Uc9Jee9F54y7eK3239uObN8P55/vRokoWIiJVpE3CALj0Un/X68qVdZzQti1ZWdnM/KKOhDFu\nHJxwAvzwh6GVsanUPhtQLAKKRUCxCF9aJYyOHf0T+WKnk6kuq0dPslatYtWWVVUPfPWVH19x993h\nFlJEJEWlTR9G1Jw5vr+6tNTPlFHDkUfy/042hl/wS84bGDNPy89/Du3bN/xRftIinHNUuArKK8rZ\n7XZTXlFOhauoskTPqbIPV+85YRwHcLga6w5X+Vmqr7fE68J8j9r+vWrdr/MbfX6WZXHHiXfUev6e\naKR3Aw0eDAccAM88Az+urV+7c2dGtOpM8eL3goSxbh089RTMnduiZU2E3RW72Va+jbJdZWzb5X/W\ntuzcvbPWZcfuHXUeiz1eXlHuv+Qrdlf5sm/IdvX1LMsi27LJycohOyu7cjvLsjAzsiyrcjGqbjfk\nnHgcNwyLTPdS27oR2a5lvdGva8I14/r+1V5XXfS8GvvjcH5sOcK4fjKdn2XJ1wCUdgkD4IYb4Npr\n/cSENWYy79yZQ3J687vFLwX7nn8eTjwR9tqrRcvZGNF5craXb2fN1jWsKVtT+XP9tvVs3LGRjds3\n+p8x65t2bGLTjk1s3bm1MhG0zm1N69zWtMptVbkeuxTkFFCQU0BeVh552VWXwtxCOhZ0JD8nv8ax\n2CU3K5fsrMiXfMyXffXt+o5Ft6v/R9KcQQHFIqBYhC8tE8bJJ/s5AF94Ac45p9rBzp3pt7sNX63/\nivXb1tOxVUefMC69NCFljbV151a+Wv8VSzYuqbF8M+sbNn+4me3l2ykqLKJrYVeKWhdRVFhEx4KO\ntM9vT1FhEft12o/2Be1pn9++8me7/HYU5hXSOrc1+dn5df7lIyJSn7Trw4h64QW45Rb45JNqE8Pe\ndhts3sxJAz/hmsOv4fR+J/sJA5cu9bPItoDyinLmrp7LjG9n8Pnqz1mwdgHz185n9dbV7NNxH/q2\n70uf9n0ql97tetOjbQ+KCoton99eX/giUi/1YTTSmWfC+PHw+utQZXLazp2htJRj+hzDe4vf4/Q1\nHf2IvxCTxfby7UxfMp1pX0/j/SXv89nKz+jdvjeH9TiMg7sdzIn7nMiALgPo16Ef2Vn1zNEuIpJA\naZswsrLg17+G3/3O3zVV+Ud5587w3Xcc0/cixr05DtZ1gaOPjvv7b9m5hakLpzL588m8U/oOBxUd\nxEn7nMTvRv2OQ3scSvuC9o26ntpnA4pFQLEIKBbhS9uEAX5Mxm23wUsv+RoHUJkwhvcczuerP2fX\n5/3IPTkOD1SKmLdmHvd/dD+TP5/MyD4j+dFBP+Lxsx/3fSUiIiksbfswol5+2c9QPnu2f2opJSVw\nySXw2Wcc9dhRvHrnt7R9agoMG9as95m7ei43vnUjHy//mCsPvZKfHfozurXp1qxriog0heaSaqLT\nTvOTzT7+eGRH27awZQsAR3T9HgXfLIUBA5p8/bJdZYx9bSzHPX4cx+99PIvHLuamUTcpWYhI2gkt\nYZjZnWY238xKzOx5M2sXc+xGM1sUOV7Pg1XjUQ6YMAFuusk/D4k2bfwkg8CR9GZ9+3xo1apJ1/5s\n5WcMe3gYa8rWsODnCxg7Yiz5OflxLH1A8+QEFIuAYhFQLMIXZg3jDeAg59wQYBFwI4CZDQQuAAYA\n3wcmWsj3iR5+OBxxBNx1F1VqGEN2dqa0XUWTrjntq2mc9MRJjD9mPE+e8ySdWnWKY4lFRJJPi/Rh\nmNnZwLnOuR+b2Q2Ac85NiBx7FbjZOfdRLa9rdh9G1JIl/kmrn8xw7L1/DuzYwe4n/sZzf7yC0z5a\nR7v8dnu+SMSLC17kin9ewZTzp3B03/jfYSUi0hyp3odxGfBKZL0nEPuc1OWRfaHq0weuuw7G/tJ8\ns9SWLWQv/5btPbrx2crPGnydD5d+yE9f+ikvX/iykoWIZJRm3VZrZtOA2N5dAxzwa+fcS5Fzfg3s\ncs493ZT3GDNmDP369QOgQ4cODBkypPJe62ibZUO3Dz20mIkTYVtuW1pt2ULxv//Nl3md2bRyFkf3\nPXqPr3/xtRf56Us/5a/X/JVDexza6PdvznZs+2xLvF8yb0f3JUt5ErldUlLC2LFjk6Y8idy+7777\nmvX9kMrbxcXFTJo0CaDy+zIMoTZJmdkY4HLgeOfcjsi+6k1SrwE3hd0kFfX667DfGQey14f/oPVN\n1/Pqcb15fv9dPHrmo3t87QXPXUDPtj2599R741qmhijWoKRKikVAsQgoFoGwmqRCSxhmdipwN3CM\nc+67mP0DgSeBw/FNUdOA/WvLDGEkDICvOx/GCyc9yLXLrmXmL85nbNnzvHfpe/W+5p1v3uGyqZcx\n/+r5FOQUxL1MIiLxkop9GA8AbYBpZjbTzCYCOOfmAc8C8/D9GleFkhXq0XtgW/49bTNbv91Ar76D\nWPjdwnrPr3AVXPvGtUw4cYKShYhkrNAShnNuf+dcX+fcsMhyVcyx251z+znnBjjn3girDHXJ7diG\na6/YwualG2jbaX+27tzKhu0b6jx/2lfTADh/4PktVcQaYtvvM51iEVAsAopF+NJ+pHet2rZlxEGb\nac9Gbn2wI/079+eL776o8/QHPn6AXwz/haYVF5GMlrEJg/XrKXDbeOzZNnSsOKDOhPH1+q/5aPlH\njB40uoULWZU68wKKRUCxCCgW4cvMhJGXB2vWYO3a8ehfjE9e70/Jstr7MabMm8J5A86jVW7Tpg8R\nEUkXmZkwsrNh3Tro0IHTToMj9j+AZ978gtq63v/5xT8544AzWr6M1ah9NqBYBBSLgGIRvsxNGN99\nB+39Q4zGX7U/a90XPPZY1dPWbVtHycoSjut3XAIKKSKSXDI7YUQey9q/a1/yuy7hhhtg7tzgtNe+\nfI1R/UYlRXOU2mcDikVAsQgoFuHL3ISxfj208xMOFrUuYvvurdx6ZxnnngubNvnT3vr6LU7Z95QE\nFlREJHlkZsLIyYFt23znN35UZK92vTj2jKUceyz85CfgHHy0/CNG9BqR4MJ6ap8NKBYBxSKgWIQv\nMxNGdjbs3Ol/RvRu35tlm5bxxz9CaSncfs9mvtnwDYO7DU5cOUVEkkjmJowdO6omjHa9WbppKQUF\nMGUK/OGJ2fQuGEhedl4CCxpQ+2xAsQgoFgHFInyZmzB27vRNUxE92vbg283fAtC3L1x83XyWfDKQ\npUvruoiISGbJ7IQRU8Moal3Emq1rglO6zmfUoIGcdRZs3ZqIQlal9tmAYhFQLAKKRfgyN2Hs2FGl\nhtG1sCtryoKEMX/tfP7rnAEMGgRjxkBF0x79LSKSNjI3YVSvYRQWVUkYX6//mv067csjj8CyZXDL\nLYkoaEDtswHFIqBYBBSL8DXrEa0pKzsbdu2q0SS1eutqAJxzLN20lN7te1OQD//4BwwfDgMHwvmJ\nm+FcRCShMreGAVWapIoKgz6MDds3kJOVQ7t8P7Cve3d48UW46iqYObPFSwuofTaWYhFQLAKKRfgy\nO2FU7/QuW4NzjmWbltGrXa8qLxk6FB56CM4+G1asaMnCiogkh8xOGDE1jFa5rcjNymXzzs2+Oapd\n7xovO+88uPxyOOssKCtrqcJ6ap8NKBYBxSKgWIQvMxNGNFHE1DAAOhR0YNOOTSzbtIyebXvW+tLx\n42H//XXnlIhknsxMGLU0SQG0yWvD5h2bWbdtHV1ad6n1pWbwl7/4O6duvjnkcsZQ+2xAsQgoFgHF\nInyZnTByqt4k1ja/LVt2bmH9tvV0bNWxzpcXFPg7p554Ap58MsyCiogkj8xOGLXVMHZuZsP2DXQo\n6FDvJbp1g6lTYexY+OCDsAoaUPtsQLEIKBYBxSJ8mZ0wqtcw8iI1jO3r6VhQdw0javBgmDQJzj3X\nz3ArIpLOMjth1NGH0ZAaRtQPfgDjxsEZZwQPXgqD2mcDikVAsQgoFuFTwojRJq9NUMOopw+jumuu\ngZEjYfRo2L07ngUVEUkemZ0wammSamgfRiwzeOABP5/h//xPPAsaUPtsQLEIKBYBxSJ8mZ0w6qph\nbGtYH0as3Fx47jl45RV4+OF4FVREJHlkdsKo5bbaxvZhxOrYEf75T/jf/4U334xHQQNqnw0oFgHF\nIqBYhC/0hGFm15lZhZl1itl3o5ktMrP5ZnZy2GWooY6R3m3y2rC6bDW52bnkZuc26dL77w/PPAMX\nXggLFza3oCIiySPUhGFmvYCTgMUx+wYAFwADgO8DE83MwixHDXU0SbXNa8u6beua/RzvUaPg9tvh\n9NPhu++adamYa46Kz4XSgGIRUCwCikX4wq5h3AtcX23fWcBk51y5c64UWAQMD7kcVdXRJNUmrw3r\nt60nN6tptYtYP/mJn6Tw3HP9s5pERFJdaAnDzM4Eljrn5lQ71BNYGrO9PLKv5dRRw8jLzqNsVxk5\nWfF5rtSECdCuHVx5JTjXvGupfTagWAQUi4BiEb5mfTOa2TSgW+wuwAHjgV/hm6OaZcyYMfTr1w+A\nDh06MGTIkMqqZ/QXpNHbhYV+e+FCKC6uPD5vxjzWzltLTr+c5l0/sv3++8VceSXccMMo7r4bDj20\nedfTtt+OSpbyJHK7pKQkqcqTyO2SkpKkKk9LbhcXFzNp0iSAyu/LMJhr7p++tV3UbBDwJlCGTyK9\n8DWJ4cBlAM65OyLnvgbc5Jz7qJbruDDKx8yZ8L3vwZQpvs0oori0mIv/cTHZWdl8c803cXu7JUvg\niCPgz3+G006L22VF0oZzftBr7FJRESzO+SV2vfp2U86LxzViP0M81838TBJNYWY45+LeNxzKM72d\nc58D3aPbZvYNMMw5t97MpgJPmtk9+Kao/YCPwyhHnepoksrNymVb+bYm3VJbnz59YPJk/wCmTz6B\n3jWfzSRSyTk/CHRPy65dwVJeHs56eXnVL/DavtTr227oPoCsLP9fMrpkZQWLWfAzusRu17Xe1GON\nuUbsLTvxXM/ObnrCCEsoCaMWDl/TwDk3z8yeBeYBu4CrwqlG1KOOTu/c7FzKdpVR1Loo7m959NF+\nZtsf/QiKi/1Av8Yojmk6y3QtHYvdu2HLFr+UldW+bN1a97Hosm1bwxNBXh7k51ddCgqqbuflwcaN\nxXTvPorcXP/rHPuzrn0FBfUfj12P/QKPfonvaV9TzsmKQ2+q/o+Er0UShnNun2rbtwO3t8R716q+\nGsaubXHr9K5u3Dh47z3/1L4JE0J5C6mmosJPCrlunb/Fed06v2zY4Pdv2gSbNwfrte0rK4M2bfxS\nWAitW+956dKl6narVsEXfvUv/upLXl7VvzbrU1zsb+MWaQktVcNILnXUMHKycnC40BJGVpZ/6NKw\nYXDssY3rz9BfTp5zMHjwKObOhZUrqy6rVtVMDOvX+y/5zp2hUye/dO4M7dv7pV07/2yTdu1qLm3b\n+p+FhfH5CzgM+r0IKBbhy8yEUcdI7+jo7qaO8m6ILl180hg9GmbP9tviOQdr1/qbBBYvrvpzyZIg\nKbRpA92711wOOsgng9jk0LFj45v/RKR2mZkw6mmSAkKrYUQde6yfOuRnP/MTFjak+SFd2medgzVr\n4Isvai5ff+2bbvr29TcK9Onj14880q/vtZevDXz4YXrEIh7S5fciHhSL8GV2wqilSQqIy0jvPfn9\n7+HQQ/0zwf/zP0N/u4TYsQPmzoXPPqu6ABxwgF/69/fJs39/2GcfX3sQkeQUyjiMeAltHMaKFdCj\nB3z4IYwYUbn7283f0vOenpyw9wm8eXGcp5utxaxZcMopvmmqe/c9n5/MnPPNRh984Jfp02HBAth3\nXzjkEL8MGQIHH+xrCSISnpQah5H06rqttoWapKKGDvVzTl17LTz1VIu8ZVytWgVvvAGvv+7v1ikv\n908ePPJIX2sYOtTfESQi6SFJ7/0IWR19GJVNUiF2elf3m9/4is60afWfV31ajERwzteGxo/3A+UP\nOABeeAGOOQbefddX3J5/Hq67zo9sDytZJEMskoViEVAswpfZNYw67pJqqRoG+Hv0H3wQrroK5sxJ\nzr/Iv/zS97VMnuwHn51/Ptx3n2/N0x1IIpkjM/swNm/2N9jPnQsDB1bu3lG+g4JbCzh/4Pk8e/6z\n8X/fepx3nm/f/9//bdG3rVN5OUydCg895GsVo0f7UeqHH97wQWUikhjqw4inJKphRN11l79r6r/+\nK7Gdwjt2wGOPwR13+FtZr7zSz8+Yn5+4MolIcsjsPoxqnd5ZloVhLdqHEbX33nDJJXDzzbUfD7t9\ntrwcJk6E/fbzzyV/9ll4/33feZ1syUJt1QHFIqBYhE81jGpys3PJscSEZfx435F8zTVw4IEt977v\nvgs//7mv2fzf/8Fhh7Xce4tI6sjMPgzn/ORAS5bUmGu8zW1tuGjwRTx8xsPxf98GuO02mD/fTx8S\nti1b4Je/9LfF3nOPb3pS/4RI6gurDyMzm6Sik9rn1KxJ5GTlJKRJKurqq+HVV/00GWH69FN/a+yu\nXb7v/7zzlCxEpH6ZmTAgmIi/mtzs3IR0eke1b+87vu+8s+r+eLbPPvUUnHoq/Pa3MGmSn5U1lait\nOqBYBBSL8GV2wqilhpGbldsic0nVZ+xY3+n87bfxv/Ztt8GvfgXvvONvkxURaajMTRjRR4tV353g\nGgZAUZH/Mn/kkWBfPGbhvOUWPwDvww9h0KBmXy5hNCNpQLEIKBbhy9yE8f77fvBeNTlZOQlPGOBH\nfv/5z76PIR7uuQf+/nd46y0/TbiISGNlbsI45JBad+dm5Sa00ztq0CA/JuLFF/12c9pnX3kF7r4b\n3nwz9WfFBbVVx1IsAopF+DI3YdQhGZqkoi6/3HdKN8eiRTBmjO8TqXYHsYhIo2TmOIx6DH14KBcO\nupDrR17fou9bmy1boFcv/6VfVNT415eX+1ljx4zxt+uKSGbQOIwWkpuVPDWMNm3g9NPhmWea9vq7\n74YOHXx/iIhIcylhVJObnRx9GFEXXQRPP9349tlly/xYjkceSb8BeWqrDigWAcUifEoY1SRTDQPg\n+OP9SOwNGxr3ut/8xg8A3HvvcMolIpkneb4Zk0Sy3FYblZ8PJ5wAGzeOavBrFizwd0YtWhReuRJJ\n99sHFIuAYhE+1TCqyc1O/Ejv6s44A156qeHn33OP77eoZZiJiEiTKWFUk2xNUgA/+AG8+mpxgwbx\nrVoFzz2X3h3daqsOKBYBxSJ8ShjVJFunN/hbart3h1mz9nzupEl+mvKm3IYrIlKfUMdhmNkvgKuA\ncuBl59whHgB3AAAMYklEQVQNkf03ApdF9l/jnHujjte3+DiMBWsX0KNtD9rlJ1d7ztVXwz77wHXX\n1X/ewQfDAw/Asce2TLlEJPmk3DgMMxsFnAEMds4NBv4Q2T8AuAAYAHwfmGiWPDd+HtjlwKRLFgBH\nHw3/+lf958yZ4++mOvrolimTiGSWMJukrgTucM6VAzjn1kb2nwVMds6VO+dKgUXA8BDLkRays4v5\n17/8wwLr8vTTMHq0fzZUOlNbdUCxCCgW4Qvzq6U/cIyZ/dvM3jGz70X29wSWxpy3PLJP6lFU5Ed+\nL1hQ9zmvvebvqBIRCUOzbgcys2lAt9hdgAPGR67d0Tk3wswOA54D9mnse4wZM4Z+/foB0KFDB4YM\nGVJ5v3X0L4pM2B41ahR7713ME0/AbbfVPP7dd7BgQTHbtgEkvrzabrntqGQpT6K2o/uSpTwtuV1c\nXMykyEyl0e/LMITW6W1mrwATnHPvRrYXASOAywGcc3dE9r8G3OSc+6iWa7R4p3cyu+MOWL3aj7Oo\nbsoUeOwxP2BPRDJbynV6Ay8AxwOYWX8gzzn3HTAV+A8zyzOzvYH9gI9DLEdaKC4uZsgQKCmp/fjb\nb/sR4Zmg+l/WmUyxCCgW4QtzhNpfgcfMbA6wA7gYwDk3z8yeBeYBu4CrVI1omEMOgdmzaz9WXOyf\nnyEiEhY9DyOFOOenK//qK+jSJdi/ebMf2LdxY62PKReRDJOKTVISZ2ZwwAGwcGHV/SUlMHiwkoWI\nhEsJI0VE22cPPLBmwpg5E4YObfkyJYraqgOKRUCxCJ8SRorp37/mtOXz58OgQYkpj4hkDvVhpJgn\nnoBXX4Wnngr2jRoF48fDiScmrFgikkTUhyEA9O0LS5ZU3bdwoe/bEBEJkxJGioi2z/bpUzVhbNrk\n75LqmUGTq6itOqBYBBSL8ClhpJiePWHlSiofplRaCv36pf+EgyKSeOrDSEG9esH06b556pVX4P77\n/cSDIiKgPgyJ0a2bn1MKYNkyn0BERMKmhJEiYttnu3aFNWv8+tKl0Lt3YsqUKGqrDigWAcUifEoY\nKaioKEgYqmGISEtRwkgRsXP+xyaM5csz6w4pqBqLTKdYBBSL8ClhpKDYhLFuHXTunNjyiEhmUMJI\nEbHts0VFQaf3unXQqVNiypQoaqsOKBYBxSJ8ShgpqHoNI9MShogkhsZhpKDp0+H66+H99yE/H3bu\n1MA9EQloHIZUKiyErVthwwZo107JQkRahr5qUkRs+2w0YWRqc5TaqgOKRUCxCJ8SRgqKJoz16zMz\nYYhIYqgPIwVt3OhHdz/7LNx3n+aREpGq1IchlQoLoazM1zA6dEh0aUQkUyhhpIjY9tmcHMjO9s/C\naNUqcWVKFLVVBxSLgGIRPiWMFNW6te/0zs9PdElEJFOoDyNF9eoFo0fDjh3+eRgiIlHqw5AqCgt9\nH4ZqGCLSUpQwUkT19tnCwsxtklJbdUCxCCgW4VPCSFHRhJGXl+iSiEimUB9Gijr1VFixAi68EMaN\nS3RpRCSZpFwfhpkdZmYfm9msyM9DY47daGaLzGy+mZ0cVhnSmWoYItLSwmySuhMY75wbCtwE3AVg\nZgOBC4ABwPeBiWYW90yYbtSHEVBbdUCxCCgW4QszYawA2kfWOwDLI+tnApOdc+XOuVJgETA8xHKk\npdat/WjvTEwYIpIYOSFe+wZgupndDRhwZGR/T+DDmPOWR/ZJPao/r7iw0P/MxCYpPbs5oFgEFIvw\nNSthmNk0oFvsLsAB44FfAL9wzr1gZucBjwEnNef9JBBNGKphiEhLaVbCcM7VmQDM7O/R4865KWb2\naOTQcqB3zKm9CJqrahgzZgz9+vUDoEOHDgwZMqTyL4lom2UmbMe2z44aNSqSMIr54guAxJevJbej\n+5KlPIncLikpYezYsUlTnkRu33fffRn9/TBp0iSAyu/LUDjnQlmAT4FjI+snADMi6wOBWUAesDfw\nJZHbe2u5hhPvnXfeqbL9wAPOgXOvvJKY8iRS9VhkMsUioFgEIt+dcf9eD20cRuQ22gcjiWE7cJVz\nblbk2I3AT4BdwDXOuTfquIYLq3yp7q9/hcsug7feguOPT3RpRCSZhDUOI7ROb+fcJ8DhdRy7Hbg9\nrPfOBJnc6S0iiaGpQVJEbPs9ZHand/VYZDLFIqBYhE8JI0VlcsIQkcTQXFIp6uOP4fDDYf58OPDA\nRJdGRJJJys0lJeFSDUNEWpoSRoqoqw8jEzu91VYdUCwCikX4lDBSlGoYItLS1IeRosrKfNLYtAna\ntk10aUQkmagPQ6po1co3R2Vik5SIJIYSRoqo3j5rBvPmZWaTlNqqA4pFQLEInxJGCtt330SXQEQy\nifowRETSjPowREQkoZQwUoTaZwOKRUCxCCgW4VPCEBGRBlEfhohImlEfhoiIJJQSRopQ+2xAsQgo\nFgHFInxKGCIi0iDqwxARSTPqwxARkYRSwkgRap8NKBYBxSKgWIRPCUNERBpEfRgiImlGfRgiIpJQ\nShgpQu2zAcUioFgEFIvwKWGIiEiDqA9DRCTNqA9DREQSqlkJw8zOM7PPzWy3mQ2rduxGM1tkZvPN\n7OSY/cPMbLaZfWFm9zXn/TOJ2mcDikVAsQgoFuFrbg1jDvBD4N3YnWY2ALgAGAB8H5hoZtHq0UPA\nT5xz/YH+ZnZKM8uQEUpKShJdhKShWAQUi4BiEb5mJQzn3ELn3CKgelvZWcBk51y5c64UWAQMN7Pu\nQFvn3IzIeX8Dzm5OGTLFhg0bEl2EpKFYBBSLgGIRvrD6MHoCS2O2l0f29QSWxexfFtknIiJJLmdP\nJ5jZNKBb7C7AAb92zr0UVsGkqtLS0kQXIWkoFgHFIqBYhC8ut9Wa2TvAdc65mZHtGwDnnJsQ2X4N\nuAlYDLzjnBsQ2f8j4Fjn3JV1XFf31IqINEEYt9XusYbRCLGFmwo8aWb34puc9gM+ds45M9toZsOB\nGcDFwP11XTCMDywiIk3T3NtqzzazpcAI4J9m9iqAc24e8CwwD3gFuCpmBN7VwF+AL4BFzrnXmlMG\nERFpGUk90ltERJJHUo70NrNTzWxBZHDfuESXJyxmVmpmn5nZLDP7OLKvo5m9YWYLzex1M2sfc37a\nDIY0s7+Y2Sozmx2zL26f3czyzGxy5DUfmlmflvt0jVNHLG4ys2VmNjOynBpzLC1jYWa9zOxtM5tr\nZnPM7L8j+zPu96KWWPwisj+xvxfOuaRa8EnsS6AvkAuUAAcmulwhfdavgY7V9k0A/l9kfRxwR2R9\nIDAL3+/ULxKjaA3xI+CwyPorwCmJ/mwN+OxHAUOA2WF8duBKYGJk/T/w44IS/rkbEYubgGtrOXdA\nusYC6A4Miay3ARYCB2bi70U9sUjo70Uy1jCG4/s2FjvndgGT8QMB05FRs5Z3FvB4ZP1xgoGNZ5JG\ngyGdc/8C1lfbHc/PHnutKcAJcf8QcVJHLKDmgFho2qDYlIiFc26lc64ksr4FmA/0IgN/L+qIRXTM\nWsJ+L5IxYVQf9JfOg/scMM3MZpjZTyP7ujnnVoH/pQG6RvZnwmDIrnH87JWvcc7tBjaYWafwih6K\nn5tZiZk9GtMMkxGxMLN++FrXv4nv/4lUjsVHkV0J+71IxoSRSUY654YBpwFXm9nR+CQSK5PvSojn\nZ0+1W7QnAvs454YAK4G743jtpI6FmbXB/8V7TeSv6zD/T6RaLBL6e5GMCWM5ENv50iuyL+0451ZE\nfq4BXsA3x60ys24Akerk6sjpy4HeMS+PxqWu/akonp+98piZZQPtnHPrwit6fDnn1rhI4zLwZ/zv\nBqR5LMwsB/8F+YRz7sXI7oz8vagtFon+vUjGhDED2M/M+ppZHvAj/EDAtGJmrSN/PWBmhcDJ+Nl/\npwJjIqddAkT/00wFfhS5s2FvgsGQK4GNZjbczAw/GPJFUoNRc8DnmMh6cz/71Mg1AM4H3g7tU8RH\nlVhEvhijzgE+j6yneyweA+Y55/4Ysy9Tfy9qxCLhvxeJvhugjjsETsXfFbAIuCHR5QnpM+6NvwNs\nFj5R3BDZ3wl4M/L53wA6xLzmRvzdD/OBk2P2fy9yjUXAHxP92Rr4+Z8CvgV2AEuAS4GO8frsQD5+\n8OgifDt4v0R/5kbG4m/A7MjvyAv4dvy0jgUwEtgd8/9iZuS7IG7/J9IgFgn9vdDAPRERaZBkbJIS\nEZEkpIQhIiINooQhIiINooQhIiINooQhIiINooQhIiINooQhIiINooQhIiIN8v8B8SSnQZTzW8UA\nAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x10761fc90>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "def plot_biquad(n, m, kk=k):\n",
    "    bb, aa = signal.zpk2tf(z[n::3], p[m::3], kk)\n",
    "    w,h = signal.freqz(bb, aa)\n",
    "    w = w * 48000 / (2 *np.pi)\n",
    "    pl.ylim((-100, 20))\n",
    "    pl.plot(w, 20 * np.log10(abs(h)));\n",
    "    return bb, aa\n",
    "b0,a0 = plot_biquad(0, 0, k**0.68)\n",
    "b1,a1 = plot_biquad(1, 1, k**0.22)\n",
    "b2,a2 = plot_biquad(2, 2, k**0.1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "-"
    }
   },
   "source": [
    "オーバーフローが発生しないよう、ゲインを適切に配分し、フィルタの順序を決定する。この例では青、緑、赤の順に適用する。ピークが生じる部分でオーバーフローさせないよう、先に減衰させるような順序としている。"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "得られた3つの伝達関数の係数を確認してみる"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 0.00959008 -0.01454516  0.00959008] [ 1.         -1.90655601  0.91159803]\n",
      "[ 0.2223618  -0.42566223  0.2223618 ] [ 1.         -1.9321054   0.95095094]\n",
      "[ 0.50490442 -0.98241691  0.50490442] [ 1.         -1.95762962  0.98620643]\n"
     ]
    }
   ],
   "source": [
    "print b0, a0\n",
    "print b1, a1\n",
    "print b2, a2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "### CMSIS-DSPLIB用に係数を変換\n",
    "\n",
    "IIRフィルタの実装は、CMSIS_DSPLIBで提供されているarm_biquad_cascade_df1_q15関数を用いる。その仕様に合わせて各係数を16ビット固定小数点の整数値に変換する。1を越えている係数があるので、16ビットに納まるよう1ビット右シフトした数値とするため32768>>1=16384を乗じ、数値を整数に丸める。"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 157. -238.  157.] [ 16384. -31237.  14936.]\n",
      "[ 3643. -6974.  3643.] [ 16384. -31656.  15580.]\n",
      "[  8272. -16096.   8272.] [ 16384. -32074.  16158.]\n"
     ]
    }
   ],
   "source": [
    "print np.rint(b0*16384),np.rint(a0*16384)\n",
    "print np.rint(b1*16384),np.rint(a1*16384)\n",
    "print np.rint(b2*16384),np.rint(a2*16384)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### 結果\n",
    "\n",
    "あとはこれをソースにコピーする。ただし、arm_biquad_cascade_df1_q15関数に使用する係数は、並べ方に注意する必要があるのと、定義の違いによりaの符号が反対であることに注意する必要がある。詳しくはリファレンス(1)を参照のこと。今回の例では以下のようになる。\n",
    "<code>\n",
    "q15_t bq_coeffs[] = {\n",
    "\t\t  515, 0,   -906,   515, 30977, -14714,\n",
    "\t\t 5171, 0, -10087,  5171, 31760, -15739,\n",
    "\t\t16384, 0, -32182, 16384, 32165, -16253,\n",
    "};\n",
    "</code>\n",
    "\n",
    "(1)http://www.keil.com/pack/doc/CMSIS/DSP/html/group___biquad_cascade_d_f1.html"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Slideshow",
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.13"
  },
  "toc": {
   "toc_cell": false,
   "toc_number_sections": true,
   "toc_threshold": 6,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
